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Abstract 

The occurrence of vortices in atomic Bose-Einstein condensates (BEC) enables a description 
of their superfluid behaviour. In this article we present a pedagogical introduction to the vortex 
physics in trapped atomic BECs. The mechanism of the vortex nucleation in an atomic BEC is 
discussed in detail. We also discuss a recently proposed approach which treats the problem of 
vortex nucleation using a one-particle Schrodinger equation with non-local and chiral boundary 
conditions. 

PACS numbers: 03.75.Kk, 03.75.Lm, 05.30.Jp, 67.40.Db 
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I. INTRODUCTION 

One of the prominent features of a superfluid is the way it behaves under rotation jl| . The 
second derivative of the free energy with respect to the rotational frequency is proportional 
to the superfluid density and vanishes for a normal fluid. Therefore a way to distinguish 
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between superfluids and normal fluids is through their response to rotation. In contrast to a 
normal fluid, which rotates like a rigid body at thermal equilibrium, the thermodynamically 
stable state of a superfluid does not rotate at low enough frequency. At higher frequencies 
the angular momentum appears as vortex filaments at which the superfluid density vanishes. 
The circulation of the velocity field flow around a closed contour which encircles the vortex, 
is quantized 

as 

. This is a consequence of the properties of the macroscopic wavefunction, 
whose phase changes by an integer multiple of 2ir around the vortex filaments. The atomic 
Bose-Einstein condensates (BEC) 0, S, 0, Q| provide a system for which this superfluid 
behaviour can be studied in the weak-coupling regime. 

In this article we present a detailed introduction to the problem of vortex nucleation 
in atomic Bose-Einstein condensates which includes the recent developments in this field 



si Is, 
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following the experimental detection of vortices in BEC 
Part of the contents overlaps with those of a recent review by Fetter and Svidzinskyjl8| 
which is a much more broad based account of vortices in BEC in general. However we 
emphasized here the problem of vortex nucleation and provide also a detailed background 
theoretical framework aiming a more general readership. 

The article is organized as follows In the first few sections we described in detail the 
vortices in trapped atomic condensates, the associated quantization of circulation, the ther- 
modynamic stability of a vortex and the evaluation of the characteristic nucleation frequency 
of a vortex based on thermodynamic arguments. In Section II, we focus on a single boson 

problem in a harmonic confinement which has been detailed in a number of references (for 

fl 

example see chapter 6, ref. We discuss the results in a static as well as in a rotating 

reference frame. The rotating trap problem is important in view of the fact that the exper- 
iments on vortices in atomic condensates are generally performed in a rotating trap. We 
show why a vortex cannot be nucleated in a mechanically stable rotating harmonic trap if 
the bosons are non-interacting whereas a finite amount of interaction makes vortex states 
energetically feasible. In Section III, we describe the Bogoliubov theory which deals with 
weakly interacting bosons at very low temperature Q (see for a review and further 
developments). The extension of this theory to the tra ppe d atomic condensate has been 



studied in a number of references 
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23 



24 



25, £61 term 



29 



30 



31 



which includes 



the cases of static as well as the rotating traps and has been reviewed in j( 
show how the Bogoliubov theory leads to the Gross-Pitaevskii (GP) equation 
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J. We 
which 
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is very successful in explaining the properties of the condensate at very low temperature 
as observed in recent experiments, both with and without vortices. We also discuss the 
meaning of the condensate wavefunction (35], 3(| and how it is related to superfluidity (see 



chapter 5 of [37J). We then explain how the vortex solutions of the GP equation lead to 
;he quantization of circulation, an idea which dates back to L. Onsager and R. P. Feynman 



flfl 



This is followed by a discussion showing how the characteristic frequency of vortex 
nucleation in a trapped condensate can be found from the solution of the GP equation 

0, 0, 0, Q- 

The main aspect of the Bogoliubov theory is the description of the 
weakly interacting bosons in terms of a condensate wavefunction and of non-interacting 
quasiparticles known as the Bogoliubov excitations. We discuss these Bogoliubov excitations 
and show how they determine the stability of a solution of the Gross-Pitaevskii equations. 



We emphasize the instability of a vortex solution in a static trap I24l. l25l. l31l . |44| and clarify 
how this instability is removed when the trap is put into rotation |3fl[ 44.l45j. We stress that 
this aspect explains why vortices are not thermodynamically stable below a characteristic 
rotational frequency of the trap 



4l| which gives 



In Section V, we describe the Thomas-Fermi (TF) approximation (40! 
analytical expressions for experimentally measurable quantities while has been very success- 
ful in explaining the behaviour of a large condensate. We also discuss the correction to the 
TF approximation required to describe the condensate profile at the boundary. We describe 
the condensate profile with and without vortex under this approximation and mention how 
an analytic expression for the vortex nucleation can be derived under this approximation 
48| . We then turn to the study of the collective excitations of the condensate within 



46 



47 



the TF approximation 



27 



471 ]. For this we introduce the hydrodynamic description of the 



condensate and show how it leads to the dispersion law of the collective excitations 
We discuss in more detail the quadrupole mode 0] since it plays an important role in the 
process of vortex nucleation. 



0fl 
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11 



12, LL 



In section VIII, we discuss experiments experiments 
and compare them with the theoretical predictions. We mention how the experimental 
observations point out the limitation of the thermodynamic critical frequency of vortex 
nucleation and points towards a more local theory of vortex nucleation which we discuss 
in the next few sections. Vortices are nucleated from the surface. So to construct such 
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a local theory one needs to know the condensate profile and its collective excitations at 
the surface of the trap. To that purpose we first account for the correction to the TF 



42 



51 



approximation required to describe the condensate profile at the boundary 
Then we discuss how this correction is included into the hydrodynamic theory to find out the 
modified dispersion relation of the surface modes |53j, |5J, |55 1 . Once the dispersion relation 
of the surface modes are identified using the Landau criterion 56[ , it is possible to determine 
the rotational frequency for which a given surface excitation becomes unstable towards the 



vortex nucleation 



32 



54 



vay agrees very well with the value experimentally observed |l5j . 
In section X, we discuss a different approach we have proposed j^l to study the problem 
of vortex nucleation where one replaces the non-linear interaction term in the GP equation 
by a set of non-local and chiral boundary conditions. The problem we study is strictly 
two-dimensional. Such boundary conditions give a natural splitting of the Hilbert space of 
the problem into bulk and edge states (which are two-dimensional analogues of the surface 
states). We show that the characteristic frequency of the vortex nucleation determined in 
this way is very close to the values experimentally observed In experiments a rotated 
trap is deformed in the plane of the rotation . If the deformation in the trap is not very 
large, it can be associated with a quadrupolar density deformation and the corresponding 
velocity field. The dynamics of this quadrupolar mode in the rotating frame determines the 
vortex nucleation in the condensate by overcomin g a surface barrier. This issue has been 



551 ] . We particularly emphasize that the characteristic frequency 



~l 



investigated in detail both experimentally |L4j, as well as theoretically 



n 



0SQ 



which we also discuss briefly. We conclude our article by summarizing these discussions. 



II. SINGLE-BOSON PROBLEM 



We start our discussion by listing the results from the solutions of the time-independent 
Scrhodinger equation of a trapped boson jl^l), namely, 



H $(r) = (-^- + V tr (r))<$>(x,y,z) = E${x,y,z) 



(2.1) 



The trap can be either static or rotating. In the later case the hamiltonian in the lab-frame 
is time- dependent. A time- independent hamiltonian (H rot ) is obtained by going into the 
co-rotating frame of the trap which is related to the static hamiltonian in the lab frame 
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(H ) by 



H ro t = Hq — n ■ L 



(2.2) 



A. Anisotropic case 



If the trap is harmonic the Schrodinger equation (|2.1|) rewrites 



^ w2 , m r 2 2 , 2 2, 2 2\ 

+ yK x +uj y y +UJ z z ) 



a>{x, yi z) = E${x,y,z) (2.3) 



The most general case is one in which all the trap frequencies are different. The problem 
is separable in Cartesian coordinates, namely, $(x, y,z) = 4>i(x)4>2(y)(p3(z). The eigenfunc- 
tions and the energy eigenvalues are respectively those of a three dimensional harmonic 
oscillator and are defined by a set of three quantum numbers n = {n x ,n y ,n z }. 

$ n = ,3/2 2 4(E a « ) e 4E ^ 2 ) FT x [K.H n (aJb~ a ) (2.4) 

V ">a- 

a 

E n = ^2(n a + ^)hw a (2.5) 

a 

Here a = x,y, z. H n (a\/ba) is a Hermite polynomial of degree n and 



where a a is the harmonic oscillator length scale corresponding to the axis a. In this most 
general case of three distinct trap frequencies, the eigenfunctions and energy levels of the 
hamiltonian H rot are not simply related to those corresponding to the non-rotating case. We 
shall not discuss this case any further. Rather in the following discussion we consider only 
some special cases of interest where such a simple relation exists. 

B. Spherical oscillator (SO) 

The hamiltonian is that of an isotropic or spherical oscillator with 
The problem is separable in spherical co-ordinates (r,6,<f)). There 

h 2 ( l d 2 L 2 l 2 2 \ , 
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The eigenfunctions and eigenvalues are described by a set of three quantum numbers n = 
(n r ,l,l z ) where n r £ N; I £ N; —l<l z < I- The eigenfunctions are written as a product 
of a Laguerre polynomial (radial part) and of a Spherical harmonic (angular part). 

£ n = ftw(2n r + l + (2.9) 
where L l n is a Laguerre polynomial and Y^i z is a spherical harmonic. These are simultane- 
ous eigenfunctions of the hamiltonian and one of the three angular momentum operators. 
Conventionally, this component is identified with the z-axis namely 

L z $ n = l g $ n (2.10) 

The energy is degenerate in l z . If the spherical trap is rotated about the z-axis, with a 
rotational frequency Vt 

H ro t = H — QL Z (2-11) 
The eigenfunctions of H rot are identical to those defined in (|2.8jl but with the eigenvalues 

E?{n) = M2nr + Kl + (2-12) 

The l z degeneracy is therefore lifted by the rotation. 



C. Two-dimensional isotropic harmonic oscillator (2DIO) 

Another case of interest corresponds to uo x = uo y = u± ^ u z . This trap geometry is 
particularly adopted in most of the experiments on atomic BEC. Using cylindrical polar 
co-ordinates (r,9,z) one finds that the problem separates into a pair of harmonic oscillators 
namely an isotropic oscillator in the x-y plane and another one dimensional oscillator along 
the z-axis. The eigenfunctions can be written as a product of the eigenfunctions of a two- 
dimensional isotropic oscillator and of a one dimensional harmonic oscillator. The energy- 
eigenvalues are the sum of the corresponding eigenvalues. We also define the aspect ratio 
Xfi which gives the relative strength of the confinement of the one dimensional harmonic 
oscillator in the ^-direction to the confinement of the isotropic oscillator, namely, 

\r=— (2.13) 

8 



Since this pair of oscillators with different dimensionalities are decoupled from each other 
we shall consider only the motion of the isotropic harmonic oscillator. In (r, 9) co-ordinates 
the hamiltonian writes as 



2m 



— ( r — ) -| £ 

dr dr r 2 



+ ]-mu 2 _r 2 (2.14) 



The eigenfunctions and energy eigenvalues can be written in terms of two quantum numbers 
{n±} = (nj z ) 



b^nUn + LV. 



1/2 



2 



$ n± = „p7p J r^e^e~+- 1 F 1 (-n, \l z \ + 1; b^r 2 ) (2.15) 



E n± = (2n + \l z \ + l)hu ± (2.16) 

where b LU± = ^g L , n G N, l z E Z and iFi(a',/3'; x) is a confluent hypergeometric function. 
Those are eigenstates of the angular momentum operator L z with eigenvalues l z . If the trap 
is rotated about the z-axis, the corresponding hamiltonian is given by 

H rot = H ± -m z (2.17) 

whose eigenfunctions are identical to those in (|2.15J) with the eigenvalues 

E r n °l = (2n +\l z \- —l z + l)Hu ± (2.18) 

LOj_ 

We shall use latter the equivalent expression of the hamiltonian in the rotating frame 

Hrot = 77^(P - mA n f + \m{ul - fi 2 )r 2 (2.19) 
2m 2 

where Aq = 17 x r and f2 = (0,0, Q). This is the Landau hamiltonian of a particle in a 
transverse magnetic field = 2f2 (written in the symmetric gauge) and in a parabolic 
confinement ^/^(a; 2 — Q 2 )r 2 . 



D. Characteristic frequency of vortex nucleation 

For a spherical oscillator, with Q < u, the ground state is (n, I, l z ) = (1, 0, 0) always 
remains the ground state. Similarly the ground state of the isotropic oscillator such that 
Q < u>± is (n,l z ) = (0,0). For a set of non-interacting rotating trapped bosons, a vortex 
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is nucleated if in the co-rotating frame a state with higher angular momentum quantum 
number becomes the ground state of the system. If this happens the corresponding rotational 
frequency is the characteristic nucleation frequency of the first vortex (fi c i). According to 
this definition, in either cases of SO and 2DIO, the the first vortex (l z = 1) is nucleated 
when the energy of the l z = 1 state in the rotating frame is equal to the ground state energy 
(l z = 0). Using (f2~T2~J) and flUBD and the relation E r ° l = E n - ttl z , we therefore get 

= N (E n (n r = U = 1.4 = 1) - S„K = M = 0,i, = 0)) = ^ 

nd = Ar (g,(n=l,l, = l)-g.(n=M, = 0)) =[jjL 

For either f2 = or f2 = a;, the confinement potential is compensated by the centrifugal 
force. The trapped condensate is therefore unstable at this point. Therefore, non- interacting 
bosons in a trap cannot have stable vortex solutions at T = 0. Since experimentally vortices 
are nucleated before this limit is reached and at very low temperatures (fcgT <C fuv) one 
concludes that a finite amount of interactions among bosons leads to the vortex nucleation. 
Thus, in order to study vortex nucleation, one needs to include interaction between bosons. 
We shall first consider the limit of a weakly interacting gas within the framework of the 
Bogoliubov j^] theory. We shall assume repulsive and local interactions among the bosons. 
The next section is devoted to the description of the Bogoliubov theory. 



III. BOGOLIUBOV DESCRIPTION OF A WEAKLY INTERACTING BOSE GAS 



The Bogoliubov description of a weakly interacting Bose gas (with and without confining 



potential) has been detailed in a large number of references amy 
We start with the following interacting hamiltonian for a iV-boson system 



21 



22 



N / ft 2 \ 
H «* = E ( -^ V ' + ) + 9 E *(r< - tv) 

i=l ^ ' i,j 



HQ- 



(3.1) 



where g 



4wh a 



m is the mass of a boson, a is the s-wave scattering length which charac- 



terizes the repulsive delta-function interaction. Also 



t t / 1 2 i 22, 2 2\ 

V tr = 2 m (^ X +UJ y y +UJ z Z ) 



(3.2) 
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In a second quantized form, the grand-canonical many-boson hamiltonian is 

drft(r) (-^-V 2 + Vtr(r)\ ftr) + % f drft \r)ft \r)ftr)ftr) - ptN (3.3) 



2m ' v ' ) r v ' 2 . 

where N = j drip^ (r)ijj(r) is the number operator and \x is the chemical potential. The 
bosons field operators satisfy 

[ftr), f (r')] = 5(r - r') ; $(r), ftr')} = [ft(r), ft (r')} = (3.4) 

In general, the field operator ip(r) can be expanded as 

ftr) = J2®Ma n (3.5) 



n 



where $ n (r) are single-particle wavefunctions and are the corresponding annihilation 
operators |6l|. They obey the commutation rules: 

[S>i] = s m,n2 , [a ni ,a n2 ] = , [a^.aj^] = . (3.6) 

Bose-Einstein condensation occurs when the number of atoms of a particular single-particle 
state becomes very large: = Nq ^> 1 and the ratio Nq/N remains finite in the thermodynamic 
limit N — > oo. In the limit N — > oo the states with jV and iVo ± 1 atoms correspond to the 
same physical configuration. Thus the operators and a can be considered as c-numbers 
since their commutator is negligible in comparison to Nq. Therefore one sets a = = v^o- 
The above definition for the Bose-Einstein condensation holds for the non-interacting gas of 
bosons. 

In presence of the interaction this definition of Bose-Einstein condensation in terms of the 
macroscopic occupation of a single boson eigenstate needs to be generalized. The generaliza- 
tion is done by requiring that one of the eigenvalues of the one body density matrix,namely 
the quantity < ft{r)ip{r) > is macroscopic We shall not discuss this issue here and 
for details see chapter 1 of 0|. Thus the field operator ftr) can be written as a sum of a c 
number ty(r) = y/N~o<$>o, plus a correction Sip(r), 

ftr) =*(r) + 6ftr) (3.7) 

Here $ is the macroscopically occupied single particle state. The classical field \l/(r) is 
defined as the condensate wavefunction. It is normalized by demanding 

dr\ty(r)\ 2 = N (3.8) 
11 



such that Nq is the total number of particles in the condensate. Plugging ()3.7|) in the 
expression for the free energy (j3.3|) and neglecting terms higher than quadratic in 8ip(r) and 
its adjoint finally yields, 

F = J drV*[H - n+ ||*(r)| 2 ]*(r) 

+ J dr (**[fT -^ + ^|*(r)| 2 ]^(» , ) + ^ + (» , )[^o-/x + ^*(r)| 2 ]*(r)) 

+ 5/*^ ,(r) ^ (r) )(-^(r)-M')(2M) (3 ' 9) 

where C = — ^V 2 + T4 r — fi + 2g\^{r)\ 2 . For later use we shall denote the 2x2 matrix 
which appears in the above expression asTis (Bogoliubov hamiltonian) . 



A. Classical approximation: Gross-Pitaevskii equation 

The first term in the expression (|3.9|) involves only the condensate wave 
and its minimization gives the stationary Gross-Pitaevskii(GP) equation 

[H + g\^(r)\ 2 ]^(r)=^(r) 



unction $(r) 
(3.10) 



This is also known as the non-linear Schrodinger equation. The time dependent GP equation 

at 



(for example see chapter 7 ^|) is obtained by replacing the fi^(r) by ih— 



ih- 



dt 



[Ho + g\*(r)\ 2 }V(r) 



(3.11) 



The consistency between two equations requires that under stationary conditions \P(r, t) 
must evolve in time as e~ % s *. The Gross-Pitaevskii energy functional which can be minimized 
to give the GP equation ([3. 10)1 is given by 

~h 2 



E GP = J dr 
and the associated Lagrangian is 
Lgp = I d r 



yy(r)\ 2 + V tr (r)\y(r)\ 2 + ^(r) 
2m 2 



(3.12) 



2 I dt dt 



(3.13) 



The quantity \l/(r) is called the condensate wavefunction since it represents the coherent 
motion of all the particles in a condensate. It is different from the first quantized many body 
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wavefunction of the iV-boson system which is a function of the co-ordinates of N bosons and 
in principle can be obtained from the diagonalization of the hamiltonian of this iV-boson 
system (|3.1|) . This wavefunction contains all correlations present in the many boson system 
whereas the condensate wavefunction only represents the coherent motion of all condensed 
particles. We shall not discuss this issue any further which is detailed in a number of 
references, for example, ( chapter 7 of (35,0] chapter 5 of 



B. Healing length of the condensate 



39(. Without confinement 



We shall now define the healing length of the condensate from the GP equation. For a 
condensate without any confinement potential, the healing length is defined as the length 
over which the wave function changes appreciably from its unperturbed value under some 
perturbation. It describes the size of the vortex core and the length over which the density 
grows from zero to its bulk value in the vicinity of a boundary 
(j3.10|) rewrites 

-— V 2 + s|*(r)| 2 tf(r) =/xtf(r) (3.14) 
2m 

The minimum energy configuration corresponds to a uniform condensate density given by 
|\[>(r)| 2 = po = ~ since the kinetic energy of such a configuration vanishes. If the density 
is locally perturbed, then far away from the perturbation it approaches this constant value 



Po 



again. The wavefunction can therefore be written as \&(r) 



ip — > 1 when r — > oo. The GP equation then becomes 

h 2 



-V 2 + |^(r)h 



JL 

Po9 



Poipir) such that 



(3.15) 



2mp g 

from which we see that the kinetic energy term is appreciable only if the wavefunction varies 
beyond the length scale defined by 

1 



h 2 



2mpo9 (87rp a)5 
£ is called the healing length or the coherence length of the condensate. 



(3.16) 



C. Vortex solution of the GP equation — Quantization of circulation 

The solutions of the GP equation f|3 . 1 Uf) can generally be written as \I/(r) = \ty(r)\e tS ( r > 
(Madelung transformation). The associated current can be defined (for details, see Will ) 
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as 

j(r) = -5- [tt(r)*Vtt(r) - (r)Vtf (r)*] = |^(r)| 2 -VS(r) = p(r)-VS(r) (3.17) 

The condensate velocity v s (r) is therefore and it is proportional to the gradient of a 
scalar function. The flow of the condensate is therefore irrotational except at the points 
where the phase is singular, namely 

Vxi) s = -VxV5(r) = fl (3.18) 
m 

This implies that the condensate can be rotated only through the formation of vortices such 
that the condensate density vanishes at the vortex cores. To illustrate this point let us 
assume that the trap is a perfect cylinder such that the wavefunction \l/ is an eigenfunction 
of the operator L z . In the cylindrical co-ordinate the phase of such an eigenfunction can 
generally be written as S(z, r, 9) = k9 such that the wavefunction is given by 

V(z,r,6) = y(z,r)e iKd (3.19) 

The single-valuedness of the wavefunction requires k to be an integer. This leads to 

v s .dr = — I VSM.dr = L—6).(drr + rd99) = — (3.20) 
m J T J r mr m 

Therefore the circulation of the local velocity v s along any closed contour T which includes 
the point where the phase is singular is quantized in units of — 

as. 



D. Leading order corrections:Bogoliubov equations 

Coming to the next order correction to the condensate wavefunction, one notes that the 
first integral in ()3.9|) is a c-number which we called E g and the second integral vanishes 
identically if \P(r) is a solution of the GP equation (|3.10|) . Since the hamiltonian (H — E g ) 
(13. 9|) is quadratic in Sip* and 5ip it can be written in a diagonalized form, namely, 

H-E g = Y,Exa\a x (3.21) 

A 

through the following Bogoliubov transformations 

5^(r) = J2(Mr)ax + vl(r)a{) (3.22) 

A 

14 



such that ma and v x satisfies the following set of the eigenvalue equations 

Cu x (r) + g^ 2 (r)v x (r) = E x u x (r) 

Cv x (r)+gV* 2 (r)u x (r) = -E x v x (r) (3.23) 

Unlike Tis in (J3.9j) which (apart from the c-number E g ) defines a set of interacting bosons, 
the hamiltonian (|3.21|) describes non-interacting quasiparticles and the state of wavefunction 
\[>(r) which satisfies the GP equation is the corresponding vacuum. The wavefunction \I/(r) 
is characterized by its winding number k and for each k have its corresponding Bogoliubov 
quasiparticles. The quasiparticles are bosons and are related to the original bosons by the 
transformation (|3.22j) . The commutators of a x , a x are same as those for the ordinary bosons 
(|3.6|) . For more details see (tJ. Since the energy of this quasiparticle is defined with respect 
to the condensate energy E g , the presence of a Bogoliubov quasiparticle with negative energy 
and the normalization (u x \u x ) — (v x \v x ) = +1 implies an energetic instability for the solution 
of the GP equation. Such instability is known as the thermodynamic instability. \l/ being a 
solution of a non-linear equation there is also another kind of instability which is known as 
the dynamical instability. This instability determines whether there exists time-dependent 
fluctuations around \I/ whose amplitude exponentially grows over time and destabilizes it. 
Relation between these two types of (in) stabilities will be discussed in section flIV Al 



E. Bogoliubov theory in a rotating frame 

Stable vortex solutions are obtained in rotating trapped condensates which consist of 
bosonic atoms (atoms with integer spin). Therefore we mention the structure of the Bogoli- 

nn 

ubov theory in a rotating frame which has been studied in |3(J, |32J . In the rotating frame 
the second quantized grand-canonical hamiltonian ()3.3j) is 

F rot = F-ft(r)(n-L)4>(r)dr (3.24) 
After the Bogoliubov decomposition we obtain 

F rot = F + J dr [tf*(r)(n-L)tf(r) + 5${y){VL ■ L)*(r) (3.25) 
+ 9*(r)(n-L)5#(r) + 5ft(T)(n-L)54>(T)]. (3.26) 
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where F is given bv ((3.9(1 . By keeping only the terms which involve only the condensate 
wavefunction, we obtain the GP equation in the rotating frame, namely 

[H + g\V(r)\ 2 - f).I]*(r) = ^(r) (3.27) 

The corresponding Bogoliubov hamiltonian in the rotating frame becomes 

fc-Q-L g^ 2 (r) \ 
H B (rot) = ( 3 - 28 ) 

y-g^(r) -[c-n-LYj 

The difference between this expression and (13.9(1 is the appearance of the term — f2 ■ L in 
the diagonal elements. For simplicity let us consider the situation when a system is rotated 
about the z-axis so that the extra term is QL Z . Under complex conjugation this operator 
changes its sign making the diagonal elements different in the 7is(rot) whereas they are 
identical in Tis apart from an overall — sign. 

F. Characteristic frequency of vortex nucleation 

The solution of the GP equation in the rotating frame is used to determine the character- 
istic frequency of vortex nucleation. The principle is the following. Let us again consider the 
case of a cylindrical symmetry ( §111 G(l . One determines the energies ((3.12(1 of a condensate 
with and without a vortex by respectively setting k equal to some non-zero integer and k = 
in the expression for \F These energies are respectively denoted as E v (k) and Eq. For a 
static trap E v (k) is always greater than Eq. In a rotating trap at a given value of Q = Q CK , 
E v (k) becomes the ground state energy of the hamiltonian in ((3.27(1 . The corresponding 
angular momentum is nhN where N is the number of atoms in the condensate. Therefore 
the characteristic rotational frequency for the nucleation of a vortex of winding number k is 
given by [4]] 

«- = ^^ < 3 ' 29 > 

The ground state energy of the hamiltonian is also the minima of the thermodynamic free 
energy. Thus this frequency is also called the thermodynamic frequency of vortex nucleation 
since it is determined by minimizing the thermodynamic free energy. Using this formula 



Dalfovo et. al. \M 



42l | evaluated numerically the critical frequency of vortex nucleation. 



Their treatment of the full GP equation ((3.10(1 includes the effect of the confinement as well 

1(3 



as that of the interaction. To understand the effect of each of these terms we start with the 
case where both the interaction and the confinement are neglected and then subsequently 
include their effect in steps. 



Non-interacting case 



To describe the non-interacting case 



we choose again the geometry of an infinitely long 
cylinder. There the superfluid flow is irrotational ()3.18|) . if the magnitude of velocity varies 

as 

v s = -9 (3.30) 
r 

where c is a constant. The above equation is valid for large distances. For small distances 
i. e. when r becomes comparable to the healing length £ (J3.16j) . v s varies rapidly and this 
definition is no longer valid. The quantization of circulation f ^lll Cj) determines the constant 
c to be 

v s = — (3.31) 

mr 

If the radius of the cylinder is b, the kinetic energy per unit length along the cylinder axis 
is given by 

T7i / \ f\ , P™ 2 2 2 f b dr irp 2 2 b 

EJk) = / 2ixrdr v. = — n n — = — k a log- (3.32) 

Jo 2 m J r m £ 

The integral is divergent. It is regulated by using a finite cut-off £ at the lower limit of r. 
This is because in presence of vortex the density p vanishes at the origin at a rate faster than 
r. The total angular momentum density per unit length of the cylinder is pnb 2 K,h. Setting 
k — 1 and using (J3.29j) one gets the thermodynamic frequency for the nucleation of the first 
vortex j3| 

• Uniform condensate (unconfined) with interaction 

Now let us consider an uniform gas of interacting bosons without confining potential Vt r . 
The geometry is again a cylinder with radius b. We evaluate the energy of the condensate 
in the presence of a vortex relative to the vortex free condensate with the same number 
of particles. The condensate wavefunction is given in (J4.15|) and because of the symmetry 
along z-axis and very low temperature it is justified to consider both the condensate and 
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the vortex state correspond to k z = 0. The extra energy per unit length of the vortex with 
K = 1 flUIg) (chapter 9 of ^). 



E v 



2nrdr 



2ixrdr 




V^(r)| 2 + g\*(r) 



\ff( r ) 



- E 



+ g\V(r) 



En 



(3.34) 



Here Eq is the energy of the bare condensate. This was evaluated first by Ginzburg and 
Pitaevskii [3] by using the solution of the GP equation (|3.1()j) and the result is (for k = 1) 



h 2 ( b 
E v = np— log 1.464- 



(3.35) 



where p corresponds to the vortex free uniform density of the condensate. In comparison to 
the expression (|3.32j) here we can see an extra factor 1.464 which accounts for the interaction 
among bosons. The corresponding characteristic frequency of the first vortex nucleation is 



n 



cl 



h , 1.4646 
log 



mb 2 



(3.36) 



The calculation is detailed in the chapter 9 of [4| and the chapter 5 of 



Q c i from the GP equation 







The method adopted in 4l|, |42j by Dalfovo et. al. to evaluate the characteristic frequency 
for the vortex nucleation for a trapped non-uniform condensate from the GP equation is very 
similar to the one described for the uniform interacting condensate. They use a geometry 
with an axial symmetry f ^ll C|) . For a given number of particles in the condensate (No) 
they solve numerically the GP equation to obtain the condensate wavefunction ty(r) in 
the presence of a vortex with winding number k and then use the solution to evaluate the 
energy functional (|3.12j) . The characteristic frequency is obtained from the relation (|3.29j) . 



In their work 



4l| ^ is plotted function of N . For N = 1 the result gives back 



Q c i = cjj_, which is the result for a set of non-interacting bosons and can be obtained 
from the solutions of the linear (single boson) Schrodinger equation (|2.21jl . For JVo > 1, 
because of the finite interaction present among bosons, the ratio is smaller than 1 
making it possible for trapped condensate to have a vortex. A numerical investigation 
of the characteristic nucleation frequency and angular momentum for a multiple vortex 
configuration was subsequently done by Butts and Rokshar 
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The reason for the reduction of Q c i in the presence of interactions can be understood as 
follows In an axially symmetric geometry and using the basis of the single particle states 
having definite angular momentum (J2.15)) (for simplicity we keep the quantum numbers n 
and n z fixed) the condensate wavefunction can be expanded as 



*(r) = 2^a>i, ( 3 - 37 ) 

=0 



The coefficient Ci z has to be fixed variationally by extremizing the following GP energy 
functional (J3.12)) in the rotating frame 



E GP ({1) = f dr 



^-\V^(r)\ 2 + V tr (r)\m(r)\ 2 + ?|^(r)| 4 - V{r)*nLMr) (3.38) 



2m 1 v yi «w> w< 2 ' 
To obtain the energy of a vortex state with winding number k this extremization has to be 
done with the constraint for the vorticity of the condensate wavefunction, namely 



J2\Qfl z = n (3.39) 



lz 

For a state with no vortex the only way to satisfy this constraint is to make C\ z = for 
l z 7^ in (|3.37|) . For k ^ 0, however it is possible to make Ci z ^ for more than one value 
of l z provided it is energetically favored. For non- interacting bosons (j2.16|) . such an option 
is energetically ruled out upto Q = u±. But for a set of bosons interacting via the contact 
interaction (J3.1|) occupying more than one angular momentum orbits reduces the interaction 
energy greatly. Therefore a vortex state becomes preferable as soon as the reduction in the 
interaction energy due to such distribution offsets the gain in the single particle energy 
(j2.16|) due to the occupation of the higher l z states. As a result a vortex can be nucleated 
for Q < u±. Later we shall show how an analytical expression for Q r -\ is obtained in the 



strongly interacting limit that agrees very well with the result of 



41 



42| for large N . 



IV. GROSS-PITAEVSKII EQUATION AND ITS LINEAR STABILITY ANALY- 
SIS 



In the section £1111 PI we have discussed the thermodynamic stability of the condensate 
wavefunction \1/ and mentioned that it is related to a dynamic stability. To check the 
dynamic stability of a solution a linear stability analysis of a steady state solution of the time- 



dependent GP equation ()3.11|) has to be performed, (7J, |22 , 



27 



471 ] . Such an analysis 
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also gives the linear response of the condensate to an external perturbation. Particularly 
the dynamical instability associated with a particular type of steady state solution leads to 
the vortex nucleation as we shall see. In the following section we describe the linear stability 
analysis. 



A. Linear response analysis of the Gross-Pitaevskii equation 

Let us consider the effect of adding a weak, sinusoidal perturbation to the trap potential. 
Following Edwards et.al. j%J one can write the time-dependent GP equation (|3.11|) under 
such weak perturbation 



ih 



2 



V(r,t) (4.1) 



where f± and ou p are respectively the amplitudes and the frequency of the sinusoidal pertur- 
bation. Since it is weak one assumes that the deviation from the solution of the unperturbed 
GP equation is also going to be small and can be written under the form 

(r, t) = exp(-ifxt) [*(r) + £*(r, t)} (4.2) 

with 5^5/ = u{r) exp(—iu> p t) + v(r)* exp(iu p t). 

After inserting Eq. ()4.2|) into Eq. (|4.1j) and keeping terms which are linear in 
u(r),v(r), f±(r) one gets 




¥(r)/ + (r) 
**(r)/!(r) 



(4.3) 



The Bogoliubov hamiltonian Tis therefore appears in the linear stability analysis. The equa- 
tions which determine the time evolution of the small fluctuations u and v are identical to 
the equations (J3.23|) (after f±(r) is set to zero and after changing p into A). The linearized 
evolution of the non-condensed part of the Bose-field operator 5ip(r) which has been dis- 
cussed earlier (section §111 Dj) is therefore formally equivalent to the linearized response of 
the condensate to a classical perturbation. Using Eq. (|3.11j) . 



(4.4) 
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The time evolution of these modes is given by e~ luJxt . This factor remains bounded in time 
provided that the imaginary part of uj\ is negative. This leads to the dynamic stability 
condition that for all uj\ Im(u\) < 0. Since the condensate wavefunction in a steady state 
is obtained by finding the local minima of the GP energy functional it can be shown that 
the dynamic stability criterion is more stringent and requires j3| Im{uo\) = 0. Using the 
normalization condition for Bogoliubov wavevectors u\ and v\ one can write the condensate 
response to an arbitrary small perturbation f(r) of frequency u p as 22 1 




9x I M r ) 
h(uj\ - Up) v v\\r) 



where 



9\ = / drf(r) [¥(r)<(r) + **(rK(r)] 



(4.5) 



(4.6) 



In the next section we shall describe how these normal modes of the Bogoliubov hamiltonian 
are calculated for a set of interacting bosons without any confinement which gives the usual 
Bogoliubov dispersion law. This is a standard text-book material and is given in many 
references. For example see 



B. Condensate in a box 



The atoms are trapped in a cubic box of size L, and we assume periodic boundary 
conditions. Solutions of the stationary Gross-Pitaevskii equation (J3.10|) can be written as 
plane waves with vanishing momentum, 

1 



$ (r 



L 3 /2 



exp(ik ■ r)|fe= 



(4.7) 



with the chemical potential \i = gN( 



Pog and the density of the condensate atoms po 



N /L 3 . In the absence of any confinement the Bogoliubov hamiltonian is translational 
invariant. Therefore one seeks its eigenvectors in the form of plane waves. 



ikr 



\ V *J 



(4.8) 
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H.B can be written in a block-diagonal form with each block (2x2 matrix) corresponding 
to a given wave vector k. 



n B [k] 



( 



2 1,2 



2m 



Po9 



-Po9 



Po9 

h 2 k 2 



2m 



Po9 



\ 
J 



(4.9) 



This matrix can be diagonalized, giving one eigenvector with the positive eigenvalue e& and 
one eigenvector with negative eigenvalue — e& where 



m i (h 2 k 2 



2po9 



1/2 



(4.10) 



2m \ 2m 

The spectrum (|4.10|) was first derived by Bogoliubov. Here only the case for the repulsive 
bosons (g > 0) is considered. For attractive bosons (g < 0) the spectrum exhibits a different 
behaviour. The Bogoliubov spectrum of repulsively interacting Bose gas strongly differs 
from those of a free Bose gas particularly for small k. Here the spectrum is linear and 
neglecting the quadratic term in ()4.10|) one gets 



e fc ~ hk 



Po9 



in 



(4.11) 



This corresponds to a spectrum of sound-waves propagating with a sound velocity c s given 
by 

(4.12) 



dk 



h dk 



Po9 
m 



It is this linear part of the spectrum which differentiates the behaviour of an interacting 
Bose gas from that of a free Bose gas. For large k the spectrum is quadratic like that of a 
free Bose gas. 



C. Thermodynamic instability of a vortex solution 

Dodd et. al. 



24j carried out the linear stability analysis of the condensate with a vortex 
(with n= 1) and find out that the corresponding Bogoliubov spectrum differs significantly 
in terms of the excitation frequency from Bogoliubov spectrum for the condensate without 
a vortex. Apart from providing the clue for the spectroscopic signature of the presence of 
a vortex their calculation also revealed the presence of a mode with negative energy and 
positive normalization in the spectrum which suggests that the vortex solution is dynamically 
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unstable. Rokshar |25( (see also |3JJ) explained this instability by pointing out that a 
Bogoliubov quasiparticle feels an effective potential which is approximately of the form (see 
the diagonal term of He) 

V eff = V tr + 2g\V(r)\ 2 -n (4.13) 

where \l/(r) is the condensate wavefunction with a vortex. Since the condensate density 
vanishes at the core of a vortex and the confinement potential also reaches it minimum 
there, V e ff can be negative in the core (relative to the chemical potential) forming a bound 
core state. The quasiparticles are also bosons. So the occupation of this bound state can 
be macroscopic and all the particles which form a vortex can be transferred to this core 
state. As a result the vortex collapse. For non-interacting bosons this can be easily verified 
since the core state is nothing but the ground state which has always a lower energy than 
the vortex state for a confined system (see the discussion in sec. flll Dl) . This has been 

n 

verified by Dodd et. a/.|24| by showing that in the limit of zero interaction, the energy 
corresponding to this mode approaches the value —hu±. The minus sign occurs since the 
energy is measured with respect to the vortex energy. Therefore in a static trap a vortex 
state is always thermodynamically as well as dynamically unstable. 

Before pointing out how this instability can be removed we briefly describe how Bogoli- 
ubov modes are calculated for a trapped condensate with or without a vortex. This has been 
done by several groups 0, 0, Q 0, Q using various schemes to study the stability 



properties of the condensate and its vortices. 

Let us assume that the atomic cloud is trapped in a cylindrical trap. This happens when 
the aspect ratio (J2.13J) \r <C 1 such that the confinement potential is given by 

V tr {r) = -mu 2 ± (r 2 + X 2 R z 2 ) « -mu 2 ± r 2 (4.14) 

and periodic boundary conditions are imposed along the z direction. The condensate wave 
function is written in the form: 

$(r) = #(r, 9, z) = m(r)e iKd e ik * z (4.15) 

The winding number k f ^IH Cj) is finite when a vortex is present and vanishes in the absence 
of a vortex. Since the trap is axisymmetric 

u x (r) = u x (r)e ik * z e i(h+K)e (4.16) 
u A (r) = v x (r)e ik * z e i( - l *^ e . (4.17) 
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The quantum numbers n = (n, l z , k z ) may take the following values: n G N; l z G Z\ and 
k z = q^-, q G Z, where L is the length of the cylinder. 

First we solve the GP equation with the ansatz 1)4.15)1 and then we diagonalize the 
Bogoliubov matrix Tis with that solution after expanding u\(r) and v\(r) in a suitable 
basis. This has been detailed in a number of references (for example see (3, 13] )• Since 
the normal modes are eigenfunctions of the angular momentum operator L z , according to 
(j4.17j) if u\(r) is an eigenfunction of L z with a particular eigenvalue l u z (in the unit of h), 
then v\(r) will be also an eigenfunction of L z with eigenvalue l u z — 2k. 

We shall now describe how we can stabilize a vortex solution. Following Fetter and 
Svidzinsky Q we denote the anomalous mode frequency as ... One way of staging 
this mode is to rotate the condensate 18, 30J. Let us consider an axisymmetric condensate 
in rotational equilibrium at an angular velocity SI around the z-axis. Since in the rotating 
frame, the hamiltonian is H — QL Z , and the Bogoliubov amplitudes have frequencies u)\(fl) = 
u)\—l z Q, where u)\ is the frequency in the non-rotating frame and l z is the angular momentum 
quantum number. It has been found in j^j that for the anomalous mode l z = —1. The 
resulting frequency in the rotating frame is 

u a (n) = u a + n, (4.18) 

Since u a is negative, the anomalous frequency in a rotating frame increases linearly towards 
zero with increasing Q and particularly u a (Q) vanishes at a characteristic rotation frequency 

Q cl = -uj a = \u a \ (4.19) 

which implies that for Q > Q c i the singly quantized vortex becomes locally stable. Linn and 
Fetter |45j carried out an explicit perturbative analysis to confirm this directly from the GP 
equation. One can also stabilize the vortex solution by increasing temperature and also by 
putting disorder which will pin the vortices j^]. The dependence of the vortex stability on 
rotation and temperature are combined to give a phase diagram for the vortices {30! . 64j . 
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V. APPROXIMATION SCHEMES TO STUDY THE BEHAVIOUR OF THE CON- 
DENSATE 



The Bogoliubov theory gives a very accurate description of the trapped BEC near the zero 
temperature, in terms of the condensate wavefunction ^f(r) and a set of collective excitations. 
It also tells whether a state represented by the condensate wavefunction is stable or not. 
Present experiments are done under conditions where further simplifications of this theory 
is possible. This enables us to interpret different experimentally measurable properties in 
terms of the known system parameters, like the number of atoms in the condensate No, the 
interaction strength (a) etc. We shall now discuss these approximation schemes in some 
detail. 

For further approximation one notes that depending on the parameters such as the num- 
ber of atoms (No) present, the trap frequency etc. the role of interaction in determining the 
properties of an atomic condensate can be different j^, 41 1. For the trapped condensate 
this was first pointed out by Baym and Pethick j^J and subsequently verified numerically 
by Dalfovo and Stringari |41|. We start with the following Gaussian variational ansatz for 
an isotropic confinement (the present derivation follows a review by Cast in 

*(') = -4p75^ /2R2 < m > 

where the spatial width R is the only variational parameter. After substitution in ()3.12|) 
this gives the following mean energy per particle 

= ggygj _ 3 tfNoaJ_ 

e ~ N ~ AmR 2 + A mU K + m R* ^ ^ 

Setting the unit of length to -4= = (h/muo) 1 / 2 and that of energy to fko one obtains: 



3 

~ 6 = l 



+ 2^ ^ 



with 



X=A/--^=. (5-4) 
V 7r y/h/mu 

The kinetic energy term scales like 1/R 2 , whereas the energy due to the confinement potential 
scales as R 2 . Since N /R 3 is the density of atoms in the condensate, parameter x 
Nqo) measures the effect of the interactions on the condensate density. The case x ^ 1 
corresponds to the weakly interacting regime, close to the ideal Bose gas limit x = 0; the 
case x ^ 1 corresponds to the strongly interacting regime. 
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A. The Thomas-Fermi approximation for an isotropic confinement 



When x ^ 1 the condensate can be well described by the Thomas-Fermi (TF) approxi- 
mation. Since b u — 1, % ^ 1 — ^o a ^ 1- If iVo is now increased, to minimize the energy 
(|5.3|) R will take a larger value. The size of the condensate R becomes larger and the kinetic 
energy term becomes less important compared to other terms since 



R >VvC 

/ h \ 2 

< 1 (5.5) 



^kin _ rnR 2 I I 



-E'harm muj 2 R 2 \mcuR 2 

This result should be contrasted against that of the noninteracting case for which the ratio is 
1. Therefore the kinetic energy term in the GP equation f)3.10j) can be neglected in this limit. 
This is known as the Thomas-Fermi (TF) approximation. We emphasize that TF regime 
is called strongly interacting, not due to a larger value of the s-wave scattering length a, 
but for the dominant role played by the interaction energy in comparison to the kinetic 
energy. The Thomas- Fermi radius R which gives the size of the system in this limit, can 
then be obtained by equating the energy due to the confinement potential and the repulsive 
interaction energy. The TF approximation leads to 

fiTF^(r) ^ V tr {r)^(r) + g\V(r )| 2 *(r ). (5.6) 

^nr)= ^ TF ~J tr{r) ) 112 (5.7) 

The chemical potential can now be expressed in terms of the other system parameters from 
the normalization condition J dr\ty(r)\ 2 = N Q where the domain of integration is limited by 
the Thomas- Fermi radius R which is then determined from the boundary condition \l/(r) = 
for r > R which yields 

Htf = ^mu 2 R 2 (5.8) 

Therefore the TF approximation determines the chemical potential fiTF and the size of 
the cloud R in terms of m, u, g, No which are all known parameters. The ground state 
condensate wavefunction ^f(r) is consequently determined as a real valued quantity which 
vanishes sharply at the boundary r = R. However (see £1111 CI ) ^(r) is generally a complex 
valued function which can be written as ^f(r) = \J p[r)e lS ^ r \ 
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B. Non-isotropic traps 



For a non-isotropic harmonic confining potential f flll A|) the boundary condition ()5.8|) 
yields the shape of the condensate as being that of an ellipsoid of radius R a along the axis 

a 

Rl= 2 -^ (5.9) 

where a = x,y, z. Using the same method as in the case of an isotropic confinement the 
condensate density and the chemical potential can be determined in terms of the system 
parameters such as N , a and u a , giving respectively 



p(r) 



\*T F (r)\ 2 = - Ivtf- V tr (r) ] 6 [pr F - V tr (r) ] = p(0) ( 1 " E ^ ) ( 1 " E ^ ) 

\a a /\a a / 



15- Noa 



(5.10) 

2/5 

^% (5.11) 



2 a 



where u> = {uj x uj y uj z ) 1 ^ , a = ^= = ^/ ^r, and the mean TF radius Rtf = ]^[(-Ra) 1//3 - in 
the TF regime, since x ^> 1, the chemical potential [j,tf satisfies 

Ht F > full, (5.12) 

At the center of the TF ellipsoid (x,y,z = 0) the density ([5.10)1 is maximum and it is given 

as 

p(0) = (5.13) 
9 

Unlike for an unbounded and uniform condensate, in a trapped non-uniform condensate 
the healing length (|3.1fj|) can only be defined in terms of a density at a reference point. 
Conventionally it is chosen as p(0), the TF density at the center of the condensate. Therefore 
we set £ = £o = [8irp(0) a]^ 1 ^ 2 . This leads to the following relation between the TF chemical 
potential ()5.13|) and the healing length 

^I = 9P(0) = »tf (5.14) 

Once the healing length is defined, the relation among different length scales in the TF 
regime can be determined from the relation 

£R TF = a 2 (5.15) 
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Since in the TF regime ^> 1, one has 

° a 1 



< 1. (5.16) 



a Rtf 

C,o (j3.16|) gives the size of a vortex core. This can be checked from the GP equation. A small 
healing length £ characterizes a small vortex core. In contrast, the healing length (and the 
vortex-core radius) in the ideal Bose gas limit are comparable with a and hence with the 
size of the condensate. 

The TF approximation (|5.5jl breaks down at the boundary of the TF region where the 
density |\I/(r)| 2 drops to zero (j5.7M5.8p . The realistic wavefunction which can be obtained 
from the solution of full Gross-Pitaevskii equation 4JJ has a rounded-off tail which vanishes 
exponentially. In the large r region (near boundary) the density is very low and therefore 
the kinetic energy cannot be any longer neglected compared to the interaction energy term 
(= p(r) 2 ) in the expression (|3.12jl . The modification due to this non-negligible kinetic energy 
will be discussed briefly in a later section (see §IXBj) . 



VI. VORTEX SOLUTION IN THE THOMAS-FERMI REGIME 



We shall now study the vortex solution under the TF expression. We consider the case 
where the trap has an axis of symmetry, namely u x = u y 
the TF radius in the x — y plane is given by R x = R y 
in this case is of the form (see section §111 Cjl 



uj_ ( §11 Cp . In such a geometry 
R i . The condensate wavefunction 



\q(r)\e iS{r) 



(6.1) 



In the cylindrical polar co-ordinate (r, 9, z) one again sets S = k9, so that 



—kB 
mr 



(6.2) 



Using the harmonic oscillator length scales a± — 1/ an< ^ a z = V v^I an d the harmonic 
oscillator energy scale 1tlj± first we write the GP equation (|3.10J) in a dimensionless form 



2m 



d , d . k 
— (r — ) + - 
dr dr r 



+ -mul{r 2 + X 2 R z 2 ) + g\V(r)\ 2 



\P(r) = /i^f(r) 



(6.3) 



Due to the presence of the centrifugal term, the solution of this equation for k ^ has 
to vanish on the z-axis. From the asymptotic behaviour of the Eq. ()6.3|) near the origin 
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it can be verified that the region over which the density changes very fast to zero is given 
by the healing length £o- Therefore in this region the kinetic energy cannot be neglected. 
Away from the origin the density profile is going to be same as the TF density without 
a vortex. The TF approximation can now be applied by neglecting the radial part of the 
kinetic energy (since it is sub-dominant relative to the centrifugal term for small r) which 
gives the following TF density profile 2^1 

Wo)(i-^-J-j)e( 1 -^--J--g («,) 

where we have used the relation (J5.14)) assuming that the correction to the /j,tf due to 
the presence of a vortex can be neglected [^j]. Thus the density profile of a vortex in TF 
approximation significantly differs from the ground state density profile only over a length 
scale £o- 



A. Critical frequency of vortex nucleation in the TF approximation 

The calculation for the critical frequency under TF approximation was carried on by 



mndh et. al. 1461, Sinha 



and Feder 



The following discussion closely follows 



The method consists in evaluating the extra energy associated with the presence of a 
vortex and then divide it by the total angular momentum, both determined within the TF 
approximation. The ratio then gives the critical frequency. Once there is a finite confinement 
along the z direction the system is not exactly like a infinite cylinder. But if the aspect ratio 
(J2.13j) \r <C 1, a cylindrical profile of the trap is still a good approximation. 

First one assumes that the TF radius satisfies the condition R± < ^. To obtain the 
energy of a vortex, Lundh et. al. proposes that for such a large condensate the domain 
of integration in the transverse plane can be divided into two parts. In the region having 
a radius r; such that £ <C r; L <C R± one can use the expression (|3.35|) after replacing b for 
Ti and for the region < r < R± one obtains the energy by integrating the kinetic energy 
term (interaction is neglected since the density is low in this region). This gives the energy 
of a vortex (k = 1) per unit length as 

h 2 ( r \ 1 r R± 

E v = 7rp(0)— log 1.464-^ + - / 27irdrmp(r)v 2 (r) 
m V £o/ 2 J r . 

, ,h 2 , /0.888i? ± \ . . 

w vrp(O)— log — - (6.5) 
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In the second integral we use the ground state TF density ()5.10|) since in this domain it is 

w/oi (E3D- 

The total angular momentum for a vortex state with k = 1 in the TF approximation is 
given by 

L z = h 



r 2 



P(0) / (1 - w )27irdr 



hp(0)nR 2 ± (6.6) 



The critical frequency is therefore 



H , /0.888.R 



"^WH^J (67) 

Energy of the condensate with a vortex is always measured relatively to the energy of 
the bare condensate. The energies due to the harmonic confinement {\muj 2 1 r 2 p) for these 
two cases differ appreciably from each other for r < £ . Hence its contribution is neglected 
in the above derivation. 

The next correction to the formula of the characteristic nucleation frequency comes by 
assuming the TF radius along the z direction is finite. As mentioned earlier the shape of the 
cloud in this case is that of an ellipsoid and as a result its TF radius in the x — y plane is z 
dependent. The energy for a vortex in this case can be obtained by integrating the energy of 
a vortex state per unit length in two-dimension ()6.5|) along the z-axis from — R z to R z after 
replacing R± by R±(z) = R±(l — j^). Such a calculation yields the following expressions: 

J-, ^Trft 2 z 2 f 0.888R ± n z 2 \ 47rp(0) h 2 f 0.671R, 

(6.8) 



h 



p(0) / ( 1 - — - — ] dz2nrdr 



„-hp(0)7rRlR z (6.9) 
15 



E v 5 ft ^ 0.671fl ± ^ 

Within this approximation it is important to note that for such an axisymmetric vortex the 
critical frequency does not depend on the TF radius along the axis of symmetry, though 
however the change in the shape of the cloud influences the absolute value of the critical 
frequency. This ends our discussion about the static properties of a TF (large) condensate. 
We shall now switch over to the description of the dynamic extension of this approximation. 
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VII. HYDRODYNAMIC THEORY OF THE CONDENSATE 



The study of a large (or equivalently strongly interacting) condensate can be extended 
to the time-dependent situations by formulating a hydrodynamic theory of the condensate. 



The hydrodynamic equations 



49( describes the time evolution of the density and the velocity 



field of the condensate. Vortex nucleation in a rotating trapped condensate is a dynamic 
process. Therefore the time-evolution of the density and the associated velocity field plays 
a pivotal role in the nucleation mechanism. For a large and strongly interacting condensate 
(TF limit) these hydrodynamic equations are identical to the Euler equations of the classical 
lydrodynamics. A linear stability analysis of the hydrodynamic equations of the condensate 
49| gives a description of the collective excitations in terms of the density fluctuations and 
the associated velocity fields. 

We shall see in the following discussion that one actually gets an analytical description 
of this collective excitations using this classical hydrodynamic theory. Since the Bogoliubov 
theory also describes the collective excitations one may ask how the collective excitations 
of the the hydrodynamic theory are related to the Bogoliubov excitations. Applying the 
TF approximation directly to the Bogoliubov equations through a large N expansion Sinha 



fsee also 



271 ]) showed that the dispersion law 49( for the collective excitations obtained 



under classical hydrodynamic approximation can be reproduced. The hydrodynamic theory 
of the condensate has been discussed in a number of references and has also been reviewed 
. We discuss it in the following section briefly. 



in 



A. Time evolution of the condensate density and velocity 



To enable the construction of a hydrodynamic theory of the condensate the wavefunction 
\I/(r) can be split into a modulus and a phase 



^(r) = yjpir)^ - * 



(7.1) 



and the condensate energy (|3.12j) rewrites as 

h 2 



E GP [p,S] = Jdr 



2m ^ ^ 2m 



V trP + f P 2 



(7.2) 
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Note that here we have explicitly divided the condensate phase by h. The corresponding 
Lagrangian is 



(7.3) 



Treating p(r,t) and S(r,t) as generalized co-ordinates, one gets respectively the following 
Euler-Lagrange equations that define the full hydrodynamic theory of the condensate 



dtS+ 2^ iVS) +Vtr + P9 = K^ff 
d tP + V- [pv s ] = 0. 



with 



v. 



3 proba 



h r _,_ T . vs 



(7.4) 
(7.5) 

(7.6) 



p 2im ' m 
is the local velocity field of the condensate. The equation (|7.5|) is the usual continuity 
equation. The time evolution of the velocity field is obtained by taking the gradient of both 
sides of the equation ()7.4IK 



mdfV 



-mVg + V tr (r ) + gp(r 



h 2 V 2 y/p 
2m yfp 



= -V(-mf s 2 + /i) 
where the chemical potential is defined as 



p = V tr + gp 



h 2 v 2 ^p 



2m ^fp 

Using the expression (J5.6j) the above equation can be rewritten as 

h 2 v 2 Vp 



(7.7) 



(7.8) 



P — Ptf 



(7.9) 



2m y/p 

For a uniform BEC (Vt r = 0) at T = 0, the chemical potential p = gpo and the pressure 
p = — |^|T = fPo- Using the same expression for pressure after replacing p by p, (17. 7|) can 
be rewritten as 



md+v 



-mw s +\4 r (r) + — 

2 p 2m y/p 



which is very similar to the Euler equation for classical hydrodynamics 



md t v s 



m(v s ■ V)t> s + Vtr(r ) H 



(7.10) 



(7.11) 
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where we use V(|?m^) = m(v s - V)v s . The only exception is the so called quantum pressure 
term, given by 

/U - /jt, TF = -- — , (7.12) 



2m yfp 

the only term in the equations where h appears. The source of this term is the kinetic energy 
of the condensate. 



B. Classical hydrodynamic approximation 

If the scale of variation of p is of the order of R, the pressure term in the equation (|7.1U|) 
is of the order whereas the quantum pressure term ()7.12|) is of the order Jf R3 . Therefore 
the quantum pressure term can be neglected as long as <C 1. This implies R 3> 

which finally gives the condition 



(7.13) 



This condition is satisfied for a large TF condensate (|5.1fij) . In this limit neglecting the quan- 
tum pressure term in ()7.1()j) we get the classical hydrodynamic equation for the condensate 

m(d t + v-V)v = -V[V tr + gp}. (7.14) 

and the chemical potential p becomes Ptf- 

Finally we write the classical hydrodynamic equations in a trap rotating with an angular 
velocity f2 

d t p + V ■[p(v s -nxr)}=0 (7.15) 



md t v s 



-mv s + V t rap{r ) H 

z p 



v.. (ft x r) 



(7.16) 



C. Collective excitations from the classical hydrodynamic theory 

Using the hydrodynamic equations ()7.5|) and (|7.14|) Stringari j49( studied the collective 
excitations of a trapped condensate. For a review see also . Let us denote the 

steady state (time-independent) density and velocity satisfying (J7.5|) and (|7.14|) by p s and 0. 
We follow the method described in the section (S 11V AJ) in order to linearize these equations 
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about the steady state solutions. This yields 

^- + V-[p s 6v s ] = (7.17) 

m^f + V[5pg} = (7.18) 

which gives 

d 2 t 5p- V- [c 2 s (r)V5p] =0 (7.19) 

while defining mc 2 (f) = p s (r)g. Therefore in a non-uniform trapped condensate the col- 
lective excitations which take the form of density fluctuations propagates as sound waves 
with a position dependent sound velocity c s (r ) unlike the collective excitations in a uniform 
condensate where such a sound velocity is constant (J4.12)) . The collective mode frequencies 
can be obtained by writing 

5p(r,t) = 5p(r)e lWdt (7.20) 

Since the steady state density is same as the TF density, i.e. p s = H2L£zJ^nll ) the eigenvalue 
equations are 

-u 2 Jp = - [VV tr (r).V5p - (p TF - V tr (r))V 2 5p] (7.21) 



The solutions of these modes have been obtained by Stringari j49j. Here we mention some 
results which will be later used to explain the process of vortex nucleation in the current 
experiments. 

D. Hydrodynamic Modes for trapped condensates 

• Spherical trap 
For the spherical trap described in (£HT|) Eq. fl7.21j) becomes 

d . uj 2 



^ r ^5p-—(R 2 -r 2 )V 2 5p 
or 2 



(7.22) 



u 2 Jp 

where R is the TF radius ()5.8j) . The eigenvalues are given by 

u 2 d = (I + 3n r + 2n r l + 2n 2 r )u) 2 (7.23) 

The normal modes of the density fluctuations are given in terms of the hypergeometric 
functions and the spherical harmonics, namely 

5p(r, t) = Cr l F(-n ri l + n r + ^,l+^; ^)Y l>h (9, 0)e"^ (7.24) 
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Details of this calculation is also available in chapter 7 of ref. 0]. Here n r gives the number of 
radial nodes and I gives the polarity of the density oscillation. This result can be compared 
with the dispersion relation for noninteracting bosons in a spherically symmetric trap (|2.9|) 
which is (after subtracting the ground state energy) 

Wci(n n /) = ui{2n r + /) (7.25) 

The modes with no radial nodes (n r = 0) are called the surface excitations for which 
(|7.23|) predicts the dispersion law 

uj d = Vl u (7.26) 

. The frequency of these modes is smaller than the harmonic oscillator result luo (|2.9|) 
except for the case of I — 1 mode (dipole mode, for details see 2(|). At very low energy 
only these modes are excited. Therefore when the condensate is rotated, initially all the 
angular momentum are carried by these modes. Their time evolution under rotation actually 
determines how angular momentum is transferred from the surface to the bulk of the systems 
and leads to the nucleation of a vortex. Specially the surface mode with I = 2 and energy 
y/2oj, known as the quadrupole mode, plays a very prominent role in the vortex nucleation 
mechanism. Experimental traps are generally non-spherical but has an axis of symmetry. 
So we identify such surface modes in this geometry in the next section. 

• Axisymmetric trap 

In an axisymmetric confinement f ^II C|) Eq. (|7.21|) writes as 

U %5p = col (r|; + Xlz^j 5 P - ^(Rl - r 2 - \ 2 R z 2 )V 2 5p (7.27) 

where \r is the aspect ratio and R x = R y = (1/\r)R z = R±. Unlike the case of a spherical 
trap here all density fluctuation modes and their dispersion cannot be written in an analytic 
form. However, using the symmetry about the z-axis some specific low lying modes can be 
studied. Because of this symmetry the normal modes of density fluctuations are eigenmodes 
of the operator L z = We use 9 as polar angle in the spherical polar co-ordinate system 
and as azimuthal angle in the cylindrical polar co-ordinate system. Also for notational 
convenience here we denote l z = ±p where p G N. One class of such solutions is of the form 

5p oc z p = (x + iyf = r p Y p ,± p (6, <f>) (7.28) 
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which have the dispersion 

u 2 d =pu 2 (7.29) 

Therefore they can be identified as the surface modes. Similarly it can be verified that there 
are solutions of the form 

5p ex z{x + iyf- 1 = r p Y p>±(p - X) (7.30) 

with the frequencies given by 

wS = (P - + ^ 2 (7.31) 

These two modes are degenerate for a spherical trap and both of them surface modes. 

To identify the quadrupole modes we note that for A# ^ 1, there is no spherical symmetry 
and as a result / is not a good quantum number. Since I determines the polarity, determining 
polarity of a general mode in this trap is difficult. But as shown in (|7.28j) and ()7.3()|) there are 
eigenmodes in such a system which can be mapped on to the eigenmodes of the condensate 
in a spherically symmetric confinement and their polarity can be subsequently identified. 
In this way one gets following two expressions for the quadrupole density modulation from 
fZZSfl and (m|) H after setting p = 2, namely 

5p oc r 2 e ±2id = (x + iy) 2 , l z = ±2 (7.32) 
5p oc rze ±ie = z(x + iy), l z = ±1 (7.33) 

with their frequencies given by ()7.29|) and (|7.31|) . Let us consider a linear combination of 
the two modes described in (J7.32|) of the form 

5p oc (x 2 - y 2 ) (7.34) 

Since this mode is a superposition of two eigenmodes with equal and opposite l z , its angular 
momentum about the axis of symmetry is zero. If such a mode becomes the steady state 
solution of the system, using (|7.5jl it can be shown that the associated quadrupolar velocity 
field v s is given as oc V(xy) Since the angular momentum about the z-axis is 0, such a field 



does not contain any vorticity. However if the trap is rotated new situation arises 
discuss this in a later section. 



We 
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E. Quadrupole modes in the presence of a vortex 



In the earlier section it has been pointed out that if the confinement is isotropic in the x-y 
plane, two quadrupole modes having angular momentum ±l z are energetically degenerate. 
Though this is true for a condensate without vortices, in the presence of a vortex this 



degeneracy is lifted 
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, . This can be explained by noting that the average velocity flow 
associated with the collective oscillation can be either parallel or opposite to the vortex 
flow, depending on the sign of the angular momentum carried by the excitation. This 
provide a chance to spectroscopically detect the presence of a vortex in the condensate. 
To that purpose Zambelli and Stringari j^J, using a sum rule based approach, derived 
analytical expression which relates the splitting between the quadrupole modes with the 
average angular momentum per particle. For the details of the sum-rule approach see 3| 
and how this is applied to study the response of the quadrupole mode is given in detail 
in 5(|. We shall just mention the expression obtained by Zambelli and Stringari which 
relates how the frequency difference between two quadrupole modes with equal and opposite 
angular momentum l z is associated with the average angular momentum per particle in the 
superfluid. Within the sum-rule approach the splitting between the two frequencies can be 
written as 

2(L Z ) 7u ±K f N ay 2/5 
w_ 2 = — -r^r = — 57=- 15 (7.35) 



where cu +2 and lu_ 2 correspond to two quadrupole modes with l z = ±2 and a± = 
This result has been subsequently used by the ENS [l^ and JILA ^ group to measure the 
angular momentum of the condensate in the presence of a vortex. 

VIII. EXPERIMENTS ON VORTICES IN BEC 

In this section we describe the experiments done to detect the vortices and to study their 
properties in the atomic BEC. We also compare the results of these experiments with the 
theoretical predictions which have been described in the earlier sections. 

The first experimental observation of vortices in a trapped condensate took place in JILA 
The exp erimental scheme is guided by a theoretical model proposed by Williams and 
Holland 66] which exploits the possibility of trapping otherwise identical (bosonic) atoms in 

37 



two different internal (hyperfine) states (chapter 3 of j^). The basic idea is as follows. 
Let us assume that identical atoms in two hyperfine states are trapped in two identical 
axially symmetric harmonic oscillator potentials, whose centers are spatially separated by 
a given distance. These two traps are rotated about the common axis of symmetry which 
passes through the center of the line joining the centers of these traps, with a frequency 
Q. Simultaneously an electromagnetic field is applied that couples the two internal atomic 
hyperfine states. The coupling is parametrized by two parameters, namely the detuning 5dt 
and the Rabi frequency TZ. The detuning 8dt gives the mismatch of the frequency of the 
coupling electromagnetic field to the frequency (energy) difference between the two hyperfine 
states. And 1Z gives the rate at which the population would oscillate between these two states 
if 5 dt = 0. 

The condensate wavefunction of this two-component BEC can be written as a two- 
component spinor. The symmetry of such wavefunction is different from that of a single- 
component BEC wavefunction [l^J. For explaining the basic idea behind the vortex nucle- 
ation in this scheme, it is enough to consider the following part of the free energy functional 
in the co-rotating frame 

F part (rot) = [drJ2 **(#o - nL g )% + - ^2*2) (8.1) 

J i=l,2 

where 1 and 2 denotes two hyperfine states. We know that ( §111 F|) the energy of a vortex 
with one unit of angular momentum (ac = 1) is shifted by Ml in the co-rotating frame 
relative to its energy in the laboratory frame. When this energy shift is compensated by 
the energy mismatch hSdt due to detuning, one gets one condensate at rest in the center of 
the atomic cloud and the other in a unit vortex state around it. This has been realized in 
JILA |8]. Subsequently the central region can be removed and one gets a vortex state in a 
single-component BEC j67|. Apart from these experiments in JILA the other experiments 
to detect vortices is analogous to the Rotating Bucket experiment for the Superfluid 4 He 
38^ . In this experiment, at very low temperature the sperfluid helium will come to the 
rest in the frame of a rotating bucket when the rotational frequency reaches a characteristic 
value. At this frequency a quantized vortex is nucleated in the superfluid. We discuss similar 
experiments done for the trapped condensate in the following section. 
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A. Rotating Bucket analogue for trapped atoms 



For the trapped atom, the role of the bucket can be played by the trapping potential. 
Therefore a similar experimental set up may be constructed if one can rotate the trap in 
which the condensate is formed. Such an experimental set-up to detect vortices was first 



realized by the ENS group 



10 



111 ] where the the trap is put into rotation by stirring it with 



a laser beam. This leads to the observation of the vortex nucleation in the rotating trap. 
Subset e m plo yi n g the s M „e teehni q ue MIT group HQ a ls „ observed the nucleation 
of a large number of vortices. Vortices nucleated in this way arrange them in the form of a 
vortex lattice 6s|. In the experiment done by the Oxford group the trap is rotated by 
using magnetic fields instead of laser stirrer. Here we describe how the trap is rotated by 
using a stirring laser beam following jl3i] . 



B. Experimental realization of a rotating trap 

The magnetic trap in which the the BEC is confined is an axisymmetric harmonic po- 
tential: 

Vtr(r) = l -mu\{x 2 + y 2 ) + -mu 2 z z 2 , (8.2) 



(TF radius in ^-direction) is much larger than its diameter (TF radius in t 

nn q 

plane) [13j, |15[. In the other limit (\r ^> 1) the condensate is like a pancake [9J. 



If the aspect ratio ()2.13jl \r <C 1 the shape of the condensate is like a cigar with its height 

;he transverse 
To create 

a rotating trap the condensate is stirred with a focused laser beam propagating along the 
symmetry axis (z-axis) of the trap and toggling back and forth very rapidly between two 
symmetric positions about the center of the trap. The electromagnetic field of the laser beam 
creates an average dipole potential in the x — y plane which gives the following additional 
effective confinement 

5V tr (r) = l -mul{e x X 2 + e Y Y 2 ) (8.3) 

where X an Y are the co-ordinates in the co-rotating frame. This potential is rotated 
through an acousto-optic deflector which rotates the co-ordinates X and Y at an angular 
frequency Q producing a rotating harmonic trap. Therefore 

X = x cos(fit) + y sin(fit) Y — —x sin(fit) + y cos(fit) (8.4) 
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Defining uj xy = ^1(1 + the total confinement potential in the lab frame is given by 

(Vtr + SVtr) (r) = l -m{uj 2 x X 2 + u 2 Y Y 2 ) + l -mu 2 z z 2 (8.5) 

1 -2/2, 2\ , \2 22 

= -muj ± (x +y) + -\ R mu ± z 

{(x 2 - y 2 ) cos(21tt) + 2xy sin(2fit)) (8.6) 

with 

u=yJ(u 2 x + u*)/2 e = (u 2 x -u J 2 Y )/(u J 2 x + uj 2 Y ) . (8.7) 

This potential is stationary in the rotating frame and oscillates periodically at frequency 
2Q in the laboratory frame ()8.6|) . The stirring frequency Q is chosen in the interval (0,u±) 
so that the condensate is stable. At the upper value of this interval, the centrifugal force 
equals the transverse restoring force of the trap. Similar technique to rotate the trap is 
employed in the MIT experiments with a varying range of intensity and frequency of the 
laser stirrer. Particularly the anisotropy e in their trap is higher in comparison to the 
other experiments. The Oxford group |l7[ puts the atom in a combination of the spherical 
quadrupole magnetic field and a rapidly rotating bias field (magnetic). Together they form a 
time-averaged orbiting potential which is anisotropic in the transverse plane and also rotates 
the trap. 

C. How to detect a vortex? 

• Prom the optical image of the condensate density 

In the rotating trap set-up the vortices will be nucleated when the stirring frequency exceeds 
a characteristic value tt c ±. Usually the stirring frequency is equal to the frequency at which 
the trap is rotated. The TF density of the condensate with a vortex ()6.4j) vanishes at the 
center of the trap (vortex core). Therefore one may detect a vortex through the imaging of 
the density profile of the condensate. To this purpose the stirring potential is switched off 
adiabatically (in a time long compared to uj^ 1 ) and the condensate density profile along the 
stirring axis is imaged to detect the presence of a vortex. The radius of the vortex core is of 
the order of the healing length (|5.14jl £o = (87rp(0)a) -1 / 2 which is too small (like £o ~ 0.2 fim 
in the experiments of the ENS group) to be observed optically in the experiments. To remove 
this difficulty the time of flight technique is used in which the trap is switched-off and the 
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. In this process the vortex core also gets 



condensate is allowed to expand for some time 
expanded. The image of the enlarged condensate is then taken to detect the presence of a 
vortex. 

• Measurement of the angular momentum 

In the presence of a single vortex with winding number k — 1 each particle in the superffuid 
carries on the average extra one unit (h) of angular momentum. Therefore the presence of a 



vortex can be verified by measuring this extra angular momentum. Both ENS 



121 and JILA 



P] did this measurement to confirm the presence of a vortex in a rotating condensate. For 
this they used the result ()7.35|) obtained by Stringari and Zambelli (s(| that in a cylindrically 
symmetric geometry in the presence of a vortex , the difference between frequencies of the two 
transverse quadrupole modes u +2 and u;_ 2 , corresponding respectively to excitations with 
angular momentum l z — 2 and l z = —2 is proportional to the average angular momentum per 
oarticle ()7.35|) . In a static condensate these two modes are degenerate. In the experiment 



12j the condensate is subjected to the dipole potential created by stirring laser beam but 
with a fixed basis (Q = in Eq. 18 .4j) such that x = X and y — Y. Therefore one just 
have an elliptically deformed trap. This non-rotating dipole potential is applied for a period 
which is smaller than the quadrupolar oscillation period. The condensate is then allowed to 
oscillate freely in the pure magnetic trap. Using time-of-flight technique its image is then 
taken. This procedure was repeated in the presence as well as in the absence of a vortex. 
From the precessional frequency of the major and the minor axis of the ellipse the difference 
u + 2 — uJ-2 is determined. Using (|7.35() . the angular momentum per atom is measured from 
it, 



D. Vortex nucleation in a deformed trap 



In a static trap the confinement potential in the x — y plane is isotropic. However 
once the rotation is switched on this is no more true since the stirring beam introduces 
an anisotropy in the confinement potential f WIIIBl E g. 18.3(1 . The anisotropy ()8.7(1 can 



be set either suddenly 



15j | or adiabaticallv jl4lll7| resulting in an elliptic deformation 



of the condensate. This has important consequences for the process of vortex nucleation. 
ENS group |lj] made a detailed study of the stationary state of the condensate when this 
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deformation is fixed at a given value and the rotational frequency is ramped up adiabatically. 



Their findings agree well with the theoretical prediction by Recati et. al. 



This stationary 



state becomes dynamically unstable beyond a characteristic rotational frequency. This leads 



to the nucleation of a vortex [59j|. On the otherhand when the deformation is switched of 
suddenly the mechanism of vortex nucleation [59!- lor! is different. We discuss this in a later 
section (flXH). 



E. Rotating the normal cloud 

In another experiment by the JILA group JjJ, instead of rotating the confining potential, 
the normal component of the atomic cloud which surrounds the superluid atomic cloud is 
rotated while the trap stays at rest. To observe vortex nucleation in this set-up the non- 
condensate atoms which is in the normal state is rotated first and then they are cooled 
below the condensation temperature. As a result the superfluid fraction grows from zero to 
a finite value. However the boundary of the cloud is always formed of the normal- component 
and rotates like a rigid body. The superfluid part in the interior rotates only by creating 
a vortex. With increasing rotation the surface of the cloud gets inflated in the plane of 
rotation because of the centrifugal force. The rotational frequency can be measured from 
the change in the size of the cloud in the plane of the rotation using the following expression 

where \° R is the aspect ratio (j2.13j) of the static trap. The number of vortices in the superfluid 
is then determined from Q. 



F. How good experiments agree with the theory? 

The comparison is done by comparing the theoretically determined ()6.10j) values of the 
characteristic frequency of the vortex nucleation to the corresponding experi- 



mental results. The critical frequency of vortex nucleation is dependent on the several pa- 
rameters such as the number of particles in the condensate, stirring mechanism, anisotropy 
in the trap potential, stirring time etc. Consequently the characteristic frequency in the 
vortex nucleation as a function of the transverse confinement frequency (oo±) can vary in the 
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experiments. This ratio is ~ 0.67 in llOj which is higher than the value calculated within 

n 

TF approximation, (« .41) |46j. Moreover ENS experiments also do not show any strong 
dependence on the number of condensate particles which is also in contradiction to the TF 
prediction. For the parameters used in the MIT experiment the theoretical (jbMUj) ro- 
tational frequency of the first vortex nucleation (« 0.08u±) is also much smaller than the 
experimentally observed value (« 0.3u±). In the JILA experiment j^, the characteristic 
frequency of vortex nucleation is given by 0.32u± < Q c i < 0.38u± which is also higher than 
the TF estimation (0.2u± < Q c i < 0.25u;_i_). This is true also for the experiment done in 



Oxford 



12J. Moreover according to the observations in ENS, MIT and Oxford the vortex 



group 

mm 



nucleation in the BEC is always preceded by a strong deformation of the surface of the 
atomic cloud. Particularly in the experiments by the ENS and the Oxford [17 [ 
this deformation has been identified as a surface quadrupole mode. Experiments 
also point out mechanism of nucleation of a vortex is dependent on whether the trap is put 
into rotation adiabatically or suddenly. These experiments motivate a more local approach 
to study the process of vortex nucleation in a trapped BEC. Corresponding developments 
are discussed in subsequent sections. 



IX. THE ROLE OF SURFACE BARRIER IN VORTEX NUCLEATION 

We have seen in the last section that the thermodynamic frequency of the vortex nu- 
cleation determined within the TF approximation ()6.1U|) is usually much smaller than the 
experimentally determined characteristic frequency of first vortex nucleation. This fact and 
the other experimentally observed features points out that the existence of the vortex state 
as a global minimum of the free energy of the condensate in the rotating frame beyond a 
certain characteristic rotation frequency is a necessary but not the sufficient criterion for 
the nucleation of a vortex. The presence of a vortex in the bulk of a condensate in a static 
axisymmetric trap is associated with an extra energy. Within the TF approximation it can 
be shown that this energy is highest when the vortex is centered at the center of the trap 
and decreases monotonically as a function of distance of the vortex center from the center of 
the trap 

0,0,0- 

Finally within TF approximation (without any boundary correction) 
this energy vanishes at the boundary of the system. The vortices are therefore nucleated 
from the surfaces. The dynamic evolution of the surface modes in a rotating trap and the 
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change in the energy landscape due to the presence of a vortex at a given distance from the 
center of the trap determines the mechanism of vortex nucleation in a trapped BEC. This 
will be discussed in the following sections. The discussion will start with the description 
of the Landau criterion of superfluidity which when applied to the surface modes in 
a trapped condensate, determines at what rotational frequency such a mode will become 
energetically unstable towards the vortex nucleation in the center of a condensate. 



A. Landau criterion of superfluidity 

Landau criterion determines the condition for the superfluidity. Let us consider a particle 
of mass m! sent into the Bose-gas with an initial velocity v p . The particle can transfer energy 
and momentum to the Bose-gas only by creating an excitation of momentum hk which causes 
viscous damping. If such an excitation is produced, with momentum hk, the momentum 
of the particle is reduced by the same amount Hk. The conservation of energy therefore 
requires 

7-1 ~^~ / 2 ~^r/ 1 1 2 4-j h^k^ . . 

Eu = —mv„ \m v n — nk\ = nk ■ v n . (9.1) 

k 2 p 2m' 1 p J p 2m' K J 

To produce a finite energy excitation, \v p \ has to satisfy the condition 



\ v p\ > 



k ■ v p 



>^>c s . (9.2) 
k 



k 

The last expression follows from Bogoliubov dispersion law (j4.10J) . This implies that a 
particle with an incoming velocity smaller than the sound velocity c s can move through the 
condensate without causing any damping. Therefore there is no channel at very low temper- 
ature through which the interacting Bose gas can dissipate. This leads to its superfluidity. 
However, for an ideal Bose gas this condition gives 

i \^E k h 2 k 

W T 2m (9 - 3) 

which makes it possible to excite an ideal Bose gas with an infinitesimal v p . Therefore it 
does not show superfluidity like a repulsively interacting Bose gas. Both systems however 
shows Bose-Einstein condensation. 

We have already mentioned that for a trapped condensate the low energy excitations are 
the surface Bogoliubov modes (J7.26|) which has finite angular momentum. Therefore the 
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channel through which a trapped condensate can dissipate is these surface excitations. Dal- 
fovo and Stringari 5^] showed how the generalization of the Landau criterion can determine 
at what rotational frequency such a surface modes can be excited to nucleate a vortex in 
the condensate. Let us denote the frequency of such a surface mode as oo s (l). Then the 
rotational frequency at which the l-th surface mode can be excited is given by 

n(0 = ^ (9.4) 

This excited surface mode overcomes a potential barrier in the surface and get in the bulk in 
the form of vortex. The minimum of the rotational frequency determined in this way ()9.4|) 
therefore gives the characteristic frequency of first vortex nucleation, namely f2 c i. 

Within epical hydrodynamical approximation CO, the onantity *fL monoton.ca.ly 
decreases with increasing / and does not give a realistic value [£4| for Q c i. This is because 
the classical hydrodynamic approximation ( §VII Bj) is no more valid in the surface region 
where the kinetic energy is appreciable and the interaction energy is relatively weak. This 
problem was tackled in various ways. Dalfovo and Stringari used 5J] a sum-rule based 
approach which fixes the upper limit of the frequency uj s (1). Khawaja et. a/.j3] modified 
the classical hydrodynamic approximation by including the surface correction to the TF 
density j^lj and obtain the dispersion law for various surface width. Their treatment was 
subsequently improved by Anglin [3] through the use of a combination of analytical and 



numerical techniques and the determined critical 



requency from the Landau criterion agrees 



well with the MIT experiment 15^ . Simula et. al. 3^] obtained the surface dispersion relation 
by numerically solving Bogoliubov equations in a rotating frame. Here we provide a brief 
discussion on the role of surface excitations in the vortex nucleation process. 

B. The order-parameter at the surface of the cloud 

To find out the dispersion relation of the surface modes one needs to to take into account 
the modification of the TF approximation ( §V Ajl near the surface of a condensate. This 
was carried out by Dalfovo et. al. 1 5 1 J . They started by identifying the characteristic length 
scale over which the condensate wavefunction at the surface goes to zero, namely the healing 
length (J3.16J) at the surface. To that purpose one starts with a spherical trap of radius R 
such that Htf = V tr (R) ()5.8|) . Then by doing a Taylor expansion of V tr (r) about the point 
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R when \R — r\ = A <C R, one obtains 



V tr {r) = Htf - \VV\\ r=R A + o(2) (9.5) 

where the derivative is taken along the direction perpendicular to the surface at r = R. 
F is the modulus of the attractive external force such that F = — VVt r evaluated at r = 
R. Therefore, near the surface the confinement potential takes the form of a linear ramp 
potential. Close to the boundary, where \r — R\ <C R, the GP equation (J3.10|) takes the form 

ft 2 d 2 

m + (r - R)Fm + g q 3 = . (9.6) 



2m dr 2 

Comparing this equation with the one obtained under the usual TF approximation (|5.6j) 
we see that the difference comes from the inclusion of the term due to the kinetic energy 
(only the most dominant contribution is taken) as well as the change in the form of the 
confinement potential. The surface healing length can now be determined as 

*.=r^r (9.7) 



V h 2 

Using this length scale as the unit of the length and scaling the wavefunction \1/ by 5 s (87iay/ 2 
the equation ()9.6|) can be written in a dimensionless form giving 

- (x + * 2 )* = . (9.8) 

Here ^ = 5 s (8na) 1 ^ 2 ^I and x = The above equation can be solved to give the surface 
wavefunction and the solutions are given in jsi|. We mention here the asymptotic form 
which gives an exponentially vanishing tail for the condensate wavefunction 

Al I 2 3/ 2 



*(x - oo) ~ ^exp {-^ S/2 ) (9.9) 

The Thomas-Fermi solution of \l/ and the solution of (|9.8|) determine the behavior of the 
wave function in two distinct regions of space: the former in the interior of the cloud, the 
latter in the boundary region. A full description of the condensate wavefunction is obtained 
by matching these two-type of solutions. The correction to the TF approximation is then 
extended to the hydrodynamic theory to obtain the dispersion relation of the surface modes. 
This will be discussed now. 
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C. Surface modes from the hydro dynamic theory 



The harmonic potential can be approximated as a linear ramp potential in the surface 
region (J9.6j) . Therefore, one can reformulate the classical hydrodynamic equation for the 



collective excitations (|7.21j) in the surface o 



the condensate by replacing the harmonic 



confinement with this linear ramp potential 53| (also see chapter 7 of Q). To study the 
surface modes one introduces a local two-dimensional co-ordinate system in the surface. The 
x axis of this local co-ordinate system is then identified with the direction of the gradient 
of this linear ramp potential. Along the other direction there is no such force and hence 
Sp in this direction can be chosen as plane wave. We use 5 S (|9.7|) as the unit of length 
and t = (^pr) 1//3 as the unit of time. Then the equilibrium density (J9.8)) has the form 
p = |\J/| 2 = —x for x < and for x > (neglecting the derivative). To obtain the surface 
dispersion relation one looks for the following type of the solutions 

5p = f(qx)e qx+iqz (9.10) 

The equation (j7.21|) gives the dispersion relation as 

uj* = 2q(l + 2n),n = 0,1,2,- •• (9.11) 



6p(x, z, t) = C n L n (-2qx)e gx+igz - iujt (9.12) 

where L n (—2qx) is the Laguerre polynomial and C n is a constant. Following [5^ we shall 
now state the conditions under which this dispersion law for the surface modes becomes 
same with the one determined f ^Vll DJ) under the global hydrodynamic approximation. For 
I much greater than n, the dispersion relation (J7.23|) becomes = u; 2 /(l+2n). The spherical 
cloud density has the I dependence of the form P/"(cos#) which has the same periodicity of 
cos(W) (or sin(/6 1 )). The wave number q on the surface of this cloud is therefore given by 

q - 2n - 1 2nR - 1 (9 13 ) 

Xq R \q R 

At the boundary the restoring force is F = muj 2 R (gradient of the harmonic potential) 
Eliminating R therefore one gets the agreement uj 2 d = u 2 l(l + 2n) = 2q(l + 2n) = uo 2 s (again 
in the dimensionless unit). For / ^> n it is therefore a good approximation to replace the 
parabolic confinement with a linear ramp potential. The reason is that the characteristic 
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penetration depth for a mode is of the order ^o±l. — (2n+i)R f rom ^ e surface. For n <C / this 
is much smaller than R and the assumption that these modes are confined in the surface 
region is then satisfied. The agreement is however true for all / when n = in the relation 

(H21. 

In the above derivation only the change in the form of the confinement potential from 
harmonic to linear ramp is taken into account. But the correction due to the surface ki- 
netic energy in Eq l9.8|) which gives a non vanishing quantum pressure ()7.12|) term in 
the hydrodynamic equation is neglected. However when the healing length 5 S (J9.7)) is not 
negligible compared to the wavelength of the surface modes this term is no more negligible 
Incorporating the corrections coming from the non- vanishing surface kinetic energy through 
a variational approach Khawaja et. al. obtained the following dispersion relation (again 
in the dimensionless unit and for lowest n) 

uj 2 s = 2g + 4g 4 [- log q + 0.15] (9.14) 

This observation agrees with the results by Fetter and Feder where they have considered 
the corrections to the Thomas-Fermi description of the condensate due to the presence of 
the boundary layer near the condensate surface. 

Anglin |55( extended this analysis of surface modes with small wave-vector q by going 
beyond Bogoliubov approximation through a combination of numerical and analytical tech- 
niques. For < q < 2 the dispersion relation (in dimensionless form) turns out to be 

uj 2 s = 2q + 1.35g 3 + 0.711g 4 (9.15) 

Then the ratio is minimized to give the characteristic rotational frequency of the first 
vortex nucleation. The critical frequency determines in this way agrees very well with the 
experimentally observed value in MIT[l5j] where a strong deformation due to the stirring 
laser beams can excite a surface modes of very large / and the deformed surface region has a 
non-negligible width. For the parameters of MIT experiment using ()9.13|) Anglin has found 
that the mode with / = 18 is unstable against the vortex formation when the rotational 
frequency Q ~ 0.3u>± in an axisymmetric trap. The agreement is however not so good for 
the other experiments Q, Q]. Particularly for the experiment in JILA p] where 

the normal cloud is rotated in a static trap, vortices are nucleated before a surface mode 
is excited. For a more detailed discussion on the conditions under which the determination 
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54| and [55|. In 



of the critical frequency based on Landau criterion works well, we refer to 
a later chapter we shall discuss a different theoretical framework to understand the vortex 
uucleatiou fom the surface which is particular* re.evaut to the esperhuents QQQ. 
Before that we shall discuss how the process of vortex nucleation can be studied by applying 
a set of non-local and chiral boundary conditions to the linear Schrodinger equation. 



r)| 2 , this equation has to be solved self- 



X. ROLE OF BOUNDARY CONDITIONS IN NUCLEATING VORTICES 

In the earlier sections we have discussed the process of vortex nucleation in the frame- 
work of repulsively interacting bosons using various approximations. In this section we shall 
provide an alternative description of this process where the effect of the interaction can 
be replaced by a set of non-local and chiral boundary conditions applied to an otherwise 
non- interacting problem. To determine the condensate wavefunction \l/(r) in a confined 
geometry like that of trapped condensate, one has to impose boundary conditions. The GP 
equation is equivalent to a single boson problem in an effective one body potential 

V E , given by 

V E = V tr + g\*(r)\ 2 (10.16) 

Because of the presence of the non-linear term g 
consistently with a suitable boundary condition Z2| . The problem becomes much simpler if 
there is a way to replace the non-linear term by a suitable choice of the boundary conditions. 
Generally such a replacement is difficult. This is because there is no known mapping available 
between the boundary conditions and the effective one body potential which they aim to 
replace ( g|\l/(r)| 2 in the present problem) except for some simple cases. Therefore the choice 
of boundary conditions depend on what type of effect generated by the effective one body 
potential the boundary condition is expected to simulate. 

In the present case we aim to study the process of vortex nucleation in a confined geom- 
etry by solving the linear Schrodinger equation with a set of such boundary conditions. For 
simplicity we consider a strictly two-dimensional problem. The proposed boundary condi- 
tions are motivated by the following consideration. From the previous discussion we see that 
the dispersion law for surface excitations determines the characteristic rotational frequency 
of vortex nucleation. For a two-dimensional problem, equivalent of the surface states are 
the edge states. The proposed boundary condition therefore should be able to isolate these 
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edge states from the bulk states. The edge states should have higher angular momentum 
relative to the bulk states. If the proposed boundary conditions can ensure this then we 
want to see whether with the increasing rotational frequency it is energetically favourable 
for a state with higher angular momentum to be transferred from the edge to the bulk. This 
process can then be identified with the nucleation of a vortex. In the following section we 
shall describe a set of boundary conditions which satisfy these conditions. Details of this 
method is given elsewhere 



A. Chiral boundary conditions (CBC) 

We consider the problem of a two-dimensional isotropic harmonic oscillator in a co- 
rotating frame rotating with angular velocity Q about the z-axis. This problem has already 
been discussed in the section flll CI There effectively the boundary conditions are put at 
infinity. Now we seek the solutions of the problem in a circular domain of radius R. The 
eigenfunctions and eigenvalues in an infinite domain are already given in (|2.15jl and (j2.18j) 
of section ^] Using (J2.19j) the current density can be derived for a given eigenstate 

fy nrn 

j = ^— (K,i. v *W. - *^V*; i: - 2i-A n |*n,d 2 ) (10.17) 

We have here changed the notation from to ^> n since they represent the condensate 
wavefunction and not a single boson wavefunction. The hamiltonian and the Schrodinger 
equation are however same as those in the single boson problem. With the type of eigenfunc- 
tions given in (j2.15|) it can be checked that the radial component of such a current density 
vanishes while its azimuthal component is given by 

je = -(--^r)\* n rf (10.18) 
m r n 

We define bn = This has the dimension of and is related to b UJ± through b n = -^b Ll}± 
Then with each angular momentum quantum number one can associate a length defined as 

a = ( 10 - 19 ) 
V o n 

For r < ri z , jg is positive and for r > r\ z it is negative while it vanishes at r = r\ z . Using 
the fact that the hamiltonian (|2.19|) is same as that of a charged particle in an effective 
magnetic field 2Qz we call the currents in these two regions respectively paramagnetic and 
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diamagnetic. For a domain of radius R we also defined quantities those are similar to the 
magnetic flux in the corresponding Landau problem. In a dimensionless form they are given 
by 

$ n = b n R 2 , <5> = b u}± R 2 (10.20) 

Now let us define the bulk and the edge regions using this particular value of ri g as 
a reference for a given angular momentum state such that the current associated to that 
particular angular momentum is respectively paramagnetic and diamagnetic in the bulk 
and the edge. Alternatively, one can define the bulk and the edge states for a disc of 
size R. The bulk states will have angular momentum l z < b^R 2 whereas the edge states 
have l z > bfiR 2 . We propose a set of non-local and chiral boundary conditions for the 
present problem which split the Hilbert space into a direct sum of two orthogonal, infinite 
dimensional spaces corresponding to bulk and edge states which respectively have positive 
and negative chirality on the boundary. The chirality is determined by the direction of the 



the boundary of the disc has eigenvalues given by 



azimuthal velocity projected on the boundary. The azimuthal velocity projected on 



X(R) = hl z - b n R 2 ) = hl z - $n) (10.21) 
The chiral boundary conditions are defined in the following way: 

1. For A > 0, namely for < &n < l z , 

dr*i,\ R = (10.22) 
This holds for any n and henceforth we shall drop the subscript n in 

2. For A < 0, namely for l z < <3>n, 

3 if) 

(fr + ri§ + bar) * l > ]r = R = (10 ' 23) 

For the first set of wavefunctions which accounts for the edge states we use Neumann bound- 
ary conditions. We could have used as well Dirichlet boundary conditions. However unlike 
Neumann boundary conditions they give an unphysical discontinuity [3, [zi | . These wave- 
functions are more and more localized towards the outer side of the system with increasing 
rotational frequency. For states with l z < $>q, whose wavefunctions are localized well inside 
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the disc, we impose the mixed boundary conditions (|10.23|) . These boundary conditions 
are akin to the boundary conditions introduced by Atiyah, Patodi and Singer (APS) in 
their studies of Index theorems for Dirac operators with boundaries Similar boundary 
conditions have also been applied to the Landau problem on manifolds with boundaries 

mm]. 



B. Spectrum with CBC and the vortex nucleation 

To describe the spectrum let us introduce the dimensionless form of the energy, namely 

e = 2mE ^ R2 = (10.24) 
nr rkj± 

where $ is defined in (J10.20)) . According to the chiral boundary conditions when $ is 
increased at a fixed — , the sign of the eigenvalues \(R) changes from positive to negative. 
Correspondingly the energy e ()10.24|) of a state with a given n and l z changes. This change 
in energy describes the corresponding transfer of a state from the edge to the bulk Hilbert 
spaces at the point $q = l z (jlO.22H10.23j) . as shown in FigO For large $, the infinite plane 
solutions (J2.16|) are reached asymptotically. 



C. Spectrum with CBC 

For an infinite system we have = n + |(l + (l — ^-)i 2 ). Therefore for a given is 
a linear function of I with slope (1 — ^-). When chiral boundary conditions are applied, this 
behaviour is approximately obeyed for the bulk states. But for the edge states the energy 
increases non-linearly with increasing angular momentum. 

The FigO shows the effect of an increase of the rotational frequency on the spectrum 
under the choice of CBC. We have plotted the energies of the bulk and the edge states 
for four different values of For each value of the quantity <3>n is increased by unit 
steps from 1 to 3 and the corresponding bulk and edge energies are shown. Under these 
conditions, the slope of the bulk energy levels increases while the slope of the edge energy 
levels goes down. The opposite behaviour is observed when, for a fixed $n, the ratio ^ 
increases. There are therefore certain values of the ratio ^-(< 1), at which the edge states 
for a given $^ = l z just intersect the bulk state corresponding to $^ = l z + 1. The first 
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£2/(0=0.75 



a3 




FIG. 1: Energy levels with chiral boundary condition. Energy levels for the first few angular 
momentum states are shown for n = and n = 1. Each curve corresponds to a given value of 
the angular momentum. They correspond to n = and l z between and 5 (thin lines for bulk 
states and thin dotted lines for edge states) as well as n = 1 and l z between — 1 and 3 with l z = — 1 
corresponds to the lowest curve (thick lines for bulk states and thick dotted lines for the edge states). 
For l z > the spectrum has a kink at the point l z = &q = — $. 

vortex is nucleated for l z — 1. The characteristic rotational frequency is denoted by Q c ±. 
The rotational frequency at which the edge states for $q = 1 intersects the bulk states of 
the i-th vortex with l z — 1 is denoted by fi C j. 

A finite Bose- Einstein condensate without any vortex is defined when the domain of radius 
R (bulk region) contains only the state l z = such that R = r\ (J10.19)) . The boundary is 
given by $q = 1 and the corresponding bulk and edge energy levels are shown as a function 
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FIG. 2: Effect of a faster rotation. In these figures we have plotted ^ as function ofl z for a set of 
— values (given above each figure). The three set of plots in each figure correspond to $n = 1, 2, 3 
(the thinnest one for = 1 and the thickest one for &q = 3). The dotted part corresponds to edge 
states while the continuous part corresponds to bulk states. For = 1 bulk states for all three 
values of <3>q fall on the same line. 

of the ratio in the Fig|21 When $q = 2 instead for $^ = 1, the state with l z — 1 is 
transferred from the edge to the bulk Hilbert space. For Q > Q cl the energy Jr(Zz) (jl0.24J) 
of any state with l z > 2 is less if <3>q = 2 instead of $q = 1. At Q = Q c i the boundary 
is given by = 2. The state with l z — 1 is now transferred to the bulk Hilbert space at 
this characteristic frequency. The bulk region has now one unit(^) of extra rotational flux. 
Thus for Q = f2 cl a vortex with k — 1 is nucleated in the bulk. At this characteristic value 
for fl the region with R = r 2 (jl0.19|) is defined as a condensate with a single vortex of / = 1. 
From these boundary conditions we can thus determine the characteristic frequency of 
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nucleation of the first vortex in terms of the trap frequency. It lies in between f2 = 0.35u;_i_ 
and Q = 0.36u;_i_. We note that this characteristic frequency of nucleation is close to the 



value Q = 0.29u± that has been observed in one of the experiments 15( where u± is the 
trap frequency in the transverse direction. 

D. Vortex nucleation at ft > Q c \ 

We have so far discussed only the the nucleation of the first vortex with k — 1 and 
we have determined the characteristic rotational frequency. With a further increase of the 
rotational frequency the boundary is successively given by $n = 3, 4, 5, • ■ ■ . Correspondingly 
more than one higher angular momentum states are transferred from the edge to the bulk 
Hilbert space at a time. The characteristic rotational frequencies are respectively given by 
0.42 < ^ < 0.43, 0.46 < ^ < 0.47, 0.48 < < 0.49. For example, if the rotational 
frequency of the trap is ramped up to Q C 2 and the system is allowed to come to equilibrium 
in the co-rotating frame, the boundary is given by <3>q = 3 instead of $q = 1. Since this is 
accompanied by the transfer of two units of the rotational flux quanta to the bulk, either 
two vortices each with n = 1 are nucleated or a single vortex with k — 2 is nucleated. The 
first of these two alternative situations implies that at Q = Q ci , i number of vortices with 
k = 1 enters. Alternatively the rotational frequency may be increased adiabatically from 
to a higher value. The first vortex will enter at Q&. The second vortex with k = 1 enters 
when the bulk states for $q = 3 just intersects edge states for $^ = 2. The corresponding 
rotational frequency is relatively higher. 

Some features of the vortex nucleation can thus be explained by replacing the non-linear 
interaction with the chiral boundary conditions. Quantitative comparison of the theory 
with experimental predictions is difficult since the experimental geometry is generally more 
complicated than a pure two-dimensional disk. Particularly we have not considered the 
deformation of the the trap from the initial circular shape in the plane of the rotation 
U y Q Q. This defecation effects the process of the vortex nucleation in these 
rotating traps. In the following section we shall briefly describe how this deformation is 
incorporated in the theoretical description of the problem of vortex nucleation. 
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XI. VORTEX NUCLEATION IN A DEFORMED TRAP 



It is mentioned in section Will Dl that to rotate the trap an anisotropy f)8.7j) is introduced 
in the confinement potential by a laser stirrer. This stirrer is then rotated to create the 
rotating harmonic trap. Subsequently, beyond a characteristic rotational frequency the 
first vortex is nucleated. In this section we discuss the role of this anisotropy in the vortex 
nucleation. To that purpose, following we start by considering the steady state solutions 
of the hydrodynamic equations (J7.15H7.iTH) in such rotating traps. 



A. Stationary states in the rotating frame 

It is already known from classical hydrodynamics that for a uniform fluid in a rota ting 
elliptical cylinder, the instantaneous induced velocity field in the lab frame is given by [3, 13 

Vd = U £2 + (11-25) 

where A and B are the major and the minor axis of the ellipse. In order to solve the 
equations (|7.15H7.l'5l one may therefore propose an ansatz velocity field of the following form 



11.26) 



Here x and y are now the co-ordinates in the co-rotating frame (identical as X and Y in 
§VIIID|) . The corresponding steady-state ellipsoidal density can be written in the usual TF 
form (EHUl) 



with 



p(r) 



m 



2„2 



11.27) 



u~ = u z x + a- 2afl 



ul = ujI + a 2 + 2afl 



(11.28) 



To see why the above form of the steady state density corresponds to a quadrupole mode 
let us rewrite it as 



p{r) 



11.29) 
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In terms of the co-ordinates in the co-rotating frame, the steady state density is therefore 
a sum of the usual TF density in a static trap (after replacing uj\ with uj\ + a 2 ) and a 
quadrupolar density. The deformation of the of the trap is given by a quantity 5 such that 

Co 2 x - Col e ~ 2o; ^ 
The continuity equation ()7.15|) requires that 



<5=^f— ^ = T ^7 (11-30) 
ooi - uof. 1 + cr 



a = -5n (11.31) 

These equations are satisfied simultaneously provided 

a 2 + a(co 2 ± - 2Q 2 ) + Vteuo\ = (11.32) 

Solution of the above equation gives the stationary solutions in the Q, a plane and has 
been analyzed in detail in 58(. We mention some special features. For a symmetric trap 
potential (e = 0) if fl < a = 0. Substituting this in (J11.29J) gives the usual TF density 
of the cloud in an axisymmetric trap (|5.10[) . However for Q < ^4=, the equation gets a 
new class of solutions which correspond to a = ±-^/2f2 2 — u\. This non-zero value of a for 
e = corresponds to a spontaneous deformation of the cloud to a quadrupolar shape. The 
physical origin of this spontaneous deformation of the trap is as follows. With the increase 
in the rotational frequency the energy of the quadrupolar state with l z = +2 monotonically 
decreases to till Q = u>±. Beyond this value of Q the energy of this mode becomes negative. 
This leads to a thermodynamic instability and the system is spontaneously deformed at this 
Q. For non-zero e, the phase diagram corresponds to two branches of solutions known 
respectively as, the normal branch and the over-critical branch. We do not detail these 
results any further and now go to the issue of vortex nucleation in such a deformed trap. 

B. Adiabatic switch on of the rotation 

The above mentioned solutions are stationary solutions of the hydrodynamic equations 
in the rotating frame. So, to produce them experimentally the rotation has to be ramped up 
adiabatically for a fixed anisotropy e. In the corresponding experiments |l3 . \v\ it is found 
that at higher rotational frequency which is generally around Q = these stationary 
solutions show dynamic instability leading to the vortex nucleation. By performing a linear 
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stability analysis of the equations of rotational hydrodynamics Sinha and Castin 59j showed 
that this is indeed the case. They found that if e is very close to 0, at Q = a large 
quadrupole oscillation makes the stationary solution dynamically unstable. For larger e 
the dynamic instability occurs over a larger range of Q around Q — Solving the time 
dependent GP equation (see also 7?j]) they also verified that this dynamic instability leads 
to the vortex nucleation. Their findings are particularly in agreement with the experimental 
results in jl^ . 



C. Sudden switch on of the rotation 



The early ENS experiments 1C| [12| and the MIT experiments 15J, [l6| come under this 
class. Here the stirring potential is rotated at a fixed frequency Q and the deformation 
(ellipticity) due to the stirrer is switched on suddenly to a given value of e. In the MIT 
experiment Q| the deformation is large and consequently it excites high polarity surface 
modes. We have already discussed in ( §IX CJ) how the characteristic nucleation frequency can 



be obtained in this case 



55| by using the Landau criterion. However in the other experiments 



the deformation is relatively less and the vortex nucleation can be understood by studying 
the quadrupole modes. Here we discuss the vortex nucleation mechanism following Kraemer 
e, aU The, 8 tud y ittvo.ves the evasion of the ene rgy of a displaced vortex UW 
in the presence of a surface quadrupolar deformation under the sudden switch on of the trap 
deformation. 

We start by evaluating the extra energy and the angular momentum associated with the 
presence of a straight vortex line with k = 1 in a perfectly cylindrical axisymmetric trap 
(axis of symmetry is again the z-axis) displaced from the center of the cylinder by a distance 
d. The corresponding expressions in TF approximation are respectively 



L z (d/R ± ) = N h 



1 - 



d 



5/2 



and 



E v (d/R ± ,fi TF ) = E v (d = 0,/xtf) 



d 



1 3/2 



11.33) 



11.34) 



where E v (d/ R±, htf) is given in (|fi.8j) . According to (jll.33|) and (jll.34j) if the vortex is 
exactly located at the boundary it neither carries extra angular momentum nor energy 
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relative to a vortex free condensate. When the trap is rotated with angular velocity f2 about 
the z-axis the energy of the system in the co-rotating frame is given by 



E v (d/R±, Q, htf) = E v (d = 0, ji TF ) 



d 



3/2 



_d_ 



1 5/2 



11.35) 



From this expression one finds that a vortex at d = is energetically favorable if the 
rotational frequency satisfies Q > ^ c i(^tf) — E v (d = 0, {j,tf)/NqH which is the thermo- 
dynamic criterion of vortex stability f ^III Fj) . However, even for Q > Q v (fi) and hence if 
E v (d/R±, f2, /i) is negative at d — 0, the energy exhibits a maximum at intermediate values 
of d that lies between and R±. Therefore if the rotation is suddenly switched on at a 
value greater than Q v , because of this energy barrier a vortex cannot be nucleated into the 
center of the trap though the vortex state is a global minima of the energy. To summarize 
the situation, the energy landscape will look like a valley representing this global minima 
surrounded by this energy barrier. As suggested by the experiments in ENS, MIT and Ox- 
ford, the disappearance of the barrier at sufficiently high angular velocity is preceded by 
a quadrupolar deformation of the atomic cloud. Therefore this energy landscape changes 
in the presence of a quadrupolar mode. In an ellipticall y d eformed trap the energy of a 
quadrupole deformation in the rotating frame is given by [s 



E Q (S, — , e, n TF ) = N Q n TF 



,2 1-eS 



■■UJ± ' 



11.36) 



In ((3(3] using this expressions the total energy of a quadrupole deformed condensate in 
presence of a displaced vortex can be evaluated. The energy landscapes are then redrawn at 
various values of the rotational frequency to understand the role of quadrupole deformation 
on the previously described energy barrier under sudden switch on of the deformation of 
the trap. At the time of the sudden switch-on of the deformation, 5 = 0. This corresponds 
to e = 0. Then the system will evolve through a set of intermediate configurations before 
attaining a stable value of the deformation parameter 5. By noting the corresponding 
changes in the energy barrier it is found 6(| that at a given value for the trap deformation 
S the a saddle is created in the energy barrier such that the vortex energy along that path 
continually decreases from the surface to the center. At this point the vortex nucleates. 
This explains the vortex nucleation mechanism under the sudden switch on of the trap 
deformation. 
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XII. CONCLUSION AND OUTLOOK 



In this article we have discussed the problem of vortex nucleation in a finite quantum 
system, namely the trapped atomic condensate. We have described how the thermodynamic 
description of the vortex nucleation needs to be modified in a finite trap geometry and the 
explanation of vortex nucleation requires a detailed analysis of the condensate at its surface 
or edge (for a two-dimensional geometry). Though the local theory takes different form 
to analyze the experiments based on different mechanism to rotate the trap, the unifying 
aspect of these approaches is the role played by surface Bogoliubov excitations. There are 
many issues related to the nucleation of vortices in BEC that have not been discussed in this 
article. They include for instance the problem of vortex nucleation in attractive condensate, 
the dynamics and evolution of the vortex lattice state and the nature of the many-vortex 
state in a rapidly rotating BEC. All these problems constitute an active domain of research 
both on experimental and theoretical fronts. Also for the topic discussed, the choice of the 
references which we have selected for our discussion is not all-inclusive due to the limitation 
of space. We take this opportunity to apologize for what is left out and to mention that 
this by no means diminishes their importance. For a more complete account we again refer 
to the review by Fetter and the references therein. 
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